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ISOTROPIC  TURBULENCE  AND  INERTIAL-RANGE  STRUCTURE 
IN  THE  ABRIDGED  LHDI  APPROXIMATION 

Robert  H.  Kraichnan 
Peterborough,  New  Hampshire 

The  abridged  LHDI  (Lagrangian-History  Direct  Interaction)  closure 
approximation  is  interpreted  physically  and  used  to  analyze  energy  transfer, 
effective  eddy  viscosities,  and  Lagrangian  spacetime  statistics  in  station¬ 
ary  and  decaying  isotropic  turbulence.  The  results  are  then  specialized 
to  the  inertial  range.  Numerical  values  are  predicted  for  the  Kolmogorov 
constant,  the  asymptotic  eddy  viscosities  due  to  inertial-range  wavenumbers, 
and  the  dimensionless  constant  in  Inoue's  formula  for  the  mean  square 
change  of  Lagrangian  velocity  with  time.  Computed  curves  are  presented 
for  the  localness  of  energy  transfer,  for  Lagrangian  spacetime  structure 
functions,  and  for  Lagrangian  spacetime  acceleration-acceleration  and 
velocity-acceleration  covariances.  Inertial-range  Eulerian  spacetime 
structure  functions  also  are  computed.  The  predicted  absolute  Kolmogorov 
spectrum  in  the  inertial  and  dissipation  ranges  is  compared  with  data  of 
Grant,  Stewart,  and  Moilliet  and  of  M.  M.  Gibson. 
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1.  INTRODUCTION 


This  paper  develops  some  analytical  and  numerical  predictions  of 
recent  closure  approximations  which  involve  tracing  the  history  of  velo¬ 
city  correlations  along  the  particle  trajectories.^-  Energy  transfer,  and 
both  Lagrangian  and  Eulerian  spacetime  correlations  are  explored  for 
stationary  and  decaying  isotropic  turbulence.  Simple  physical  interpre¬ 
tations  of  the  closure  formulas  and  their  consequences  are  attempted. 

The  results  are  finally  specialized  to  the  inertial  range,  and  numerical 
values  are  computed  for  all  the  dimensionless  constants  and  universal 
functions  which  appear  in  the  formulas.  The  quantities  evaluated  are 
listed  in  the  abstract.  Analytical  and  numerical  predictions  for  two- 
particle  dispersion  are  developed  in  a  paper  to  follow. 

The  closure  approximations  of  Ref.  1  are  based  on  a  generalized 
velocity  u(x,t|r),  defined  as  the  velocity  measured  at  time  r  in  that 
fluid  element  which  arrives  at  x,  in  laboratory  coordinates,  at  time  t. 
Thus, 

u(x,t|t)  =  u(x,t),  (1.1) 

where  u(x,t)  is  the  Eulerian  velocity,  and  u(a,tQ|t)  =  w(a,t),  where  tQ 
is  the  initial  time  and  w(a,t)  is  the  velocity  of  the  particle  initially 
at  a_  (Lagrangian  velocity). 

The  full  function  u(x,t|r)  can  be  computed  in  terms  of  the  Eulerian 
velocity  by  the  equation^- 

[3/3t  +  u(x,t)*7]u(x,t|r)  =  0,  (1.2) 

with  the  boundary  condition  (1.1).  The  interpretation  of  (1.2)  is  simple. 
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The  vector  u(x,t|r)  is  the  velocity  of  a  particular  fluid  element  at  time 
r,  and  (x,t)  is  a  point  on  the  spacetime  trajectory  of  this  element.  Thus 
t  may  be  called  the  labeling  time  and  r,  the  measuring  time.  The  same 
fluid  element  can  be  labeled  by  any  other  point  (xf ,t’>  on  its  trajectory, 
and  the  velocity  measured  at  time  r  is  independent  of  which  point  is  used 
for  labeling.  Equation  (1.2)  expresses  this  by  stating  that  the  substantial 
derivative  of  u(x,t|r)  in  the  spacetime  (x,t)  vanishes  everywhere.  In 
other  words,  u(x,t|r)  at  fixed  r  is  constant  along  every  particle  trajec¬ 
tory  in  (x,t)  spacetime. 

Equation  (1.2)  thus  has  a  peculiar  character.  It  expresses  only  a 
coordinate  transformation  through  which  the  velocity  field  at  time  r  is 
labeled  by  the  laboratory  coordinates  the  fluid  elements  had  (or  will  have) 
at  time  t.  The  dynamics  are  fully  determined  by  the  equation  of  motion 
for  the  Eulerian wlocity  u(x,t).  The  role  of  (1.2)  is  to  provide  a 
redundant,  but  useful,  infinite  family  of  transformed  representations  of 
the  Eulerian  field.  Such  a  family  is  useful  because  approximate  statistical 
equations  are  the  ultimate  goal.  Although  the  Lagrangian  velocity  is  fully 
determined  by  the  Eulerian  velocity  in  any  realization,  simple  statistics 
of  the  Lagrangian  velocity  are  not  determined  by  simple  statistics  of  the 
Eulerian  velocity. 

No  restriction  to  incompressibility  is  im ^plied  by  (1.2),  which  is 
valid  whatever  the  equation  of  motion  for  u(x,t).  The  latter  will  be 
the  incompressible  Navier-Stokes  equation  in  this  paper. 

Because  the  labeling  transformation  as  well  as  the  dynamics  determine 
u(x,t]r),  the  generalized  velocity  has  properties  that  are  unexpected  from 
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experience  with  Eulerian  analysis.  It  was  pointed  out  in  Ref.  1  thai. 
divergence less  Eulerian  velocity  does  not  imply  divergenceless  u(x,t|r) 
for  r  1  t.  To  explore  another  property,  consider  a  hypothetical  flow 
maintained  in  a  statistically  steady  and  isotropic  state  at  very  low 
Reynolds  number  by  means  of  suitably  random  stirring  devices.  If  D  is  a 
characteristic  length  scale  of  the  stirring  process,  the  characteristic 
scale  of  spatial  variation  of  the  Eulerian  field  u(x,t)  will  also  be  D, 
and  the  condition  of  very  low  Reynolds  number  implies  that  variations  at 
scales  small  compared  to  D  will  be  strongly  damped.  If  a  Fourier  analysis 
of  u(x,t)  is  made,  wavenumbers  >>D-1  should  be  very  weakly  excited. 

Now  consider  u(x,t|r)  for  large  values  of  |t-r| .  Although  both 
u(x,t)  and  u(x,r)  display  spatial  variation  at  scale  D,  the  x  variation 
of  u(x,t|r)  becomes  rapid  without  limit  as  |t-r|  increases.  This  is 
because  fluid  elements  which  are  close  together  at  the  labeling  time  t 
will  typically  be  farther  apart  at  time  r.  No  matter  how  small  |x-x' |  is 
taken,  u(x,t|r)  and  u(x' ,t|r)  typically  are  velocities  measured  at  widely 
separated  parts  of  the  fluid  if  | t-r |  is  large  enough.  For  sufficiently 
large  | t-r | ,  the  Fourier  analysis  of  u(x,t|r)  shows  excitation  extending 
to  indefinitely  high  wavenumbers.  In  particular,  the  excitation  does  not 
fail  off  at  the  viscous  cutoff  wavenumber  for  the  Eulerian  velocity. 

The  assumption  of  low  Reynolds  number  puts  these  observations  in  their 
most  striking  form  but  is  otherwise  unnecessary.  To  examine  a  further 


property,  retain  the  stationarity,  isotropy,  and  homogeneity  conditions 
but  relax  the  restriction  to  small  Reynolds  number.  The  Eulerian  covariance 


(u(x,t).u(x'  ,r)>,  Where 


denotes  ensemble  average,  is  an  even  function 


V'i 

& 


of  t-r  and  has  zero  slope  at  r  =  t.  This  is  not  true  of  the  Lagrangian 
covariance  (u(x,t |t)«u(x'  ,t|r)^,  unless  x'  =  x.  The  behavior  of  the 
Lagrangian  covariance  can  be  understood  as  follows.  Two  fluid  elements 
which  are  a  distance  |x-x' |  apart  at  time  t  were  typically  a  greater 
distance  apart  at  any  sufficiently  earlier  time  r.  As  the  elements  draw 
together  under  the  straining  action  of  the  turbulence,  the  joint  effects 
of  molecular  and  eddy  viscosity  on  the  steepening  velocity  gradient 
between  the  elements  tends  to  reduce  their  velocity  difference,  while  the 
stirring  action  tends  typically  to  increase  the  difference.  If  |x-x' |  is 
intermediate  between  the  scale  characteristic  of  the  stirring  and  the  scale 
characteristic  of  dissipation,  the  effect  of  the  stirring  on  relative 
velocity  is  suppressed  compared  to  the  effect  of  molecular  and  eddy 
viscosity.  The  typical  net  result  is  a  reduction  of  the  velocity  difference. 
As  the  elements  again  separate,  their  velocities  tend  to  be  bett^er 
correlated  than  before  the  close  approach.  This  suggests  that  the  Lagrang¬ 
ian  covariance  is  an  asymmetrical  function  of  r  -  t ,  for  x  i  ,  and  can 
have  a  positive  slope  at  r  =  t.  However,  if  x  =  x* ,  the  covariance  measures 
the  autocorrelation  of  the  velocity  of  a  single  fluid  element,  whence 
stationarity  and  the  equivalence  of  all  elements  require  that  the  covariance 
be  an  even  function  of  r-t. 

These  properties  are  demonstrated  analytically  later  in  the  paper. 

The  slope  at  r  -  t  is  shown  to  be  simply  related  to  the  energy-transfer 
function.  In  the  limit  of  very  low  Reynolds  number,  transfer  is  negligible, 
stirring  and  viscosity  act  on  the  same  spatial  scale,  and  the  slope  tends 


to  zero. 
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2.  THE  CLOSURE  EQUATIONS 

In  isotropic,  homogeneous  turbulence,  the  covariance  tensor 

U^(x,t|r;x'  ,t'  |r')  =  (u.(x,t |r)u^(x' ,t'  |r'))  (2.1) 

is  fully  determined  by  the  scalars 

US(k;t|r;t'  |r' )  =  (2ir)~3fe~i~'(s~*-'  )P.j(k)Uij(x,t|r;x'  ,t'  |r*  )d(x-x* ) , 


.  UC(k;t|rjt»|r')  =  ( 2ir)“3/e~i~* )nij(k)Uij(x,t |r;x' ,t* |r' )d(x-x* ) , 

(2.2) 


where 

V*)  »  ^kj/k2,  Pij(k)  =  6£j  -  k.k^k2. 

s  c 

U  and  U  are  the  scalars  of  the  solenoidal  and  antisolenoidal  parts 
of  U..,  respectively.  U^(k;t |r;t' |r' )  vanishes  if  either  t  =  r  or  t'  =  r'. 

1 J 

The  closure  approximations  of  Ref.  1  involve  both  the  covariance 
tensor  and  a  generalized  infinitesimal  Green's  tensor  G.  _j(x,t|r;x'  ,t'  | r* ) , 
which  gives  the  average  change  in  u(x,t|r)  resulting  from  an  applied 
infinitesimal  perturbation  in  u(x',t'|r').  The  definition  and  properties 
of  the  Green's  tensor  are  discussed  in  detail  in  Ref.  1.  The  isotropic 
Green's  tensor  is  fully  determined  by  the  scalars 


GS(k;t  | r ;t *  | r* )  =  j  K^^-^^P.^^G.^x.tlrjx'.t' |r')d(x-x'), 
GC(k;t|rjt' |r»)  =  ^"^’^"^^^(kjG.^x.tlrjx' ,t'  |r')d(x-x'). 


(2.3) 


GS(kjt|rjt' | r' )  and  GC(k;t |r;t' |r' )  may  be  interpreted  as  average  infini¬ 
tesimal  response  functions  for  Fourier  mode  )c.  By  definition. 
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GS(k;t|r;t|r)  =  GC(k;t  |r;t  |r)  =  1, 

GS(k;t|r;t' |r')  =  GC(k;t|r;t* |r' )  =0  (r  <  r').  (2.4) 

Two  closure  approximations  were  proposed  in  Ref.  1.  The  first,  called 
the  LHDI  (Lagrangian-History  Direct-Interaction)  approximation,  involves 
all  the  time-argument  values  of  the  covariance  and  Green's  functions.  A 
closed  subset  of  the  final  equations  involves  only  functions  for  which 
t  =  t' .  This  approximation  makes  formidable  storage  and  time  demands  in 
numerical  computations.  The  second,  simpler  closure,  called  the  abridged 
LHDI  approximation,  involves  only  the  solenoidal,  pure  Lagrangian  functions 
U(k;t|r)  and  G(k;t|r)  defined  by 

U(k;t|r)  =  US(k;t(t;t|r)  (t  >  r) 

*  U(k;r|t)  (t  <  r), 

G(k;t|r)  =  GS(k;t|t;t|r).  (2.5) 

The  LHDI  approximation  was  obtained  in  two  stages.  First,  the  direct- 
interaction  approximations  for  triple  moments  were  constructed.  Then,  the 
formulas  were  altered  in  a  systematic  way  so  that  the  final  equations 
exhibited  the  following  properties  of  the  exact  turbulence  dynamics: 

Energy  conservation;  inviscid  equipartition  equilibrium  and  fluctuation- 
relaxation  relations;  invariance  of  Lagrangian  correlations  under  random 
Galilean  transformation;  and  invariance  of  certain  averages  under  the 
transformation  from  Eulerian  to  Lagrangian  coordinates.  The  crucial 
proj§rty  in  this  list  is  invariance  to  random  Galilean  transformation. 

This  invariance,  not  retained  by  the  unaltered  direct-interaction  approxi¬ 
mation,  is  a  formal  statement  of  the  obvious  fact  that  random,  uniform 
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translations  of  the  individual  flows  in  a  homogeneous  ensemble  cannot 
affect  the  mean  energy-transfer.  The  LHDI  approximation  appears  to  be 
the  unique  alteration,  within  the  general  direct-interaction  framework, 
which  restores  this  invariance  for  the  full  Lagrangian  covariance 
tensor  without  sacrificing  the  other  properties  listed.  No  cutoffs, 
parameters,  or  arbitrary  functions  are  invoked.  The  motivations  of  the 
LHDI  approximation,  and  the  details  of  the  alterations,  are  described  at 
length  in  Secs.  1  and  6  of  Ref.  1. 

The  abridged  closure  was  obtained  in  Ref.  1  as  a  further  approximation 
to  the  LHDI  equations.  It  appears  to  be  the  unique  alteration  of  the 
direct-interaction  equations  which  gives  a  closed  set  of  equations  for 
the  restricted  quantities  (2.5)  and  exhibits  the  consistency  properties 
listed  above.  Thus  the  abridged  LHDI  equations  may  be  considered  a 
closure  in  their  own  right,  rather  than  an  approximation  to  an  approximation. 

Both  closures  assume  that  the  velocity  field  is  multivariate  normal 
at  an  initial  instant  tQ.  In  the  limit  of  very  small  Reynolds  number, 
both  closures  yield  triple-moment  formulas  which  are  indistinguishable 
from  the  first  term  in  an  expansion  of  the  triple  moments  in  powers  of 
turbulent  Reynolds  number.1 

The  analysis  in  the  present  paper  is  based  on  the  abridged  LHDI 
equations,  with  the  exception  of  the  discussion  of  Eulerian  structure- 
functions  in  Sec.  7.  Quantitative  estimates  of  how  the  numerical 
predictions  change  in  the  unabridged  approximation  have  not  been  attempted. 
Qualitative  changes  in  the  predictions  are  not  anticipated. 

In  the  case  of  isotropic  turbulence,  the  abridged  equations  of  Ref.  1, 


Sec.  10  may  be  written 


(3/3t  +  2vk2)U(k;t|t)  = 
t 

2/J^dpdqB^J  [G(k;t |s)U(p;t|s)  -  G(p;t |s)U(k;t |s)]U(q;t | s)ds,  (2.6) 
t0 

(3/3t  t  vk2)U(k;t|r)  = 

t  t 

-U(k;t|r)JJAdpdqCkp<iJ  U(q}t|s)ds  +  JJ^dpdqD^tKpjt  |r)  J  U(q;t|s)ds 

r 

+  JJ^dpdqJ  CBjCpqG(k ;rjs)U(p;t |s)  -  D^pqG(p;r|s)U(k;t|s)]U(q;t|s)ds 
t0 
t 

-  JJAdPd<lJ  Bkpq[G(p»t |s)U(k;r|s)  -  DkpqG(k;t  |e)U(psr |s)]U(q;t |s)ds 

t0 

t 

-  J/^dpdqJ  [Bkpq  -  D^pq]G(p;t |s)U(k;t |s)U(qir|s)ds  (t  Z  r),  (2.7) 

t0 

(3/3t  +  vk2)G(k;t|r)  = 

t  t 

-  G(k;t|r)//AdpdqCkpq/  U(q;t|s)ds  +  /JAdpdqDkpqG(p;t |r)/  U(q;t|s)ds 

t 

-  //AdPdqCDkpq  -  Bkpq]G(p;t |r)U(qjt |r)/  G(k;s|r)ds 

t 

-  JJ^dpdqJ  [BkpqG(p;t|^G(k;s|r)  -  D^pqG(k;t |s)G(p;s |r)]U(q;t |s)ds 

(t  z  r),  G(k;r|r)  =  1.  (2.8) 

Here  v  is  kinematic  viscosity,  the  integration  J  j  ^  is  over  the  part  of  the 
(p,q)  plane  such  that  k,  p,  q  can  form  the  sides  of  a  triangle,  and 
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Bkpq  =  "P!q<xy  *  z3>,  Ckpq  =  .kpq(l  -  y2), 

VMp/k)2,W>  <2-9) 

where  x,  y,  z  are  the  cosines  of  the  interior  angles  opposite  the  triangle 

2 

sides  k,  p,  q,  respectively.  These  equations,  with  (2.5),  determine 
U(k;t|r)  and  G(k;t|r)  for  all  t  i  r  i  tQ  if  the  initial  values  U(k;tQ|t0) 
are  prescribed.  The  interpretation  of  (2.6)-(2.8)  will  be  developed  in 
the  Sections  to  follow. 


3.  ENERGY  TRANSFER 

The  kinetic-energy  spectrum  function 

E(k,t)  *  2nk2U(k;t|t)  (3.1) 

obeys  the  balance  equation 

(3/3t  +  2vk2)E(k,t)  =  T(k,t),  (3.2) 

2 

where  T(k,t)/(2Trk  )  is  the  right-hand  side  of  (2.6).  Conservation  of 
energy, 

00 

/  T(k,t)dk  =  0,  (3.3) 

0 


follows  from  the  identity 


2  2 
k  B,  =  p  B  ,  . 
kpq  K  pkq 


T(k,t)  can  be  written  in  the  form 
.  «  p+k 

T(k,t)=y/dpJ  dqT(k,p,q,t), 
2  0  |p-kj 


(3.4) 


(3.5) 
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2 

where  T(k,p,q,t)/(4irk  )  is  the  symmetrical  part  of  the  integrand  of  the 
p,q  integration  in  (2.6).  [The  A  integration,  although  written  asymmetrically 
in  (3.5),  is  symmetric  in  p  and  q  so  that  only  the  symmetrical  part  of  the 
integrand  can  contribute.]  The  detailed  conservation  property  of  the 
Navier-Stokes  equation, 

T(k,p,q,t)  +  T(p,q,k,t)  +  T(q,k,p,t)  =  0,  (3.6) 

follows  from  (3.4). 

The  total  rate  of  energy  transfer  into  all  wavenumbers  >k  is 

o 

00  00  p'+k* 

n(k,t)  =  i  /  dk'J  dp'/  dq'T(k',p',q',t).  (3.7) 

1  k  0  Ip'-k'l 

Since  T(k' ,p' ,q' ,t)  is  symmetrical  in  p'  and  q',  the  total  value  of  (3.7) 
must  equal  twice  the  contribution  from  q'  >  p'.  But  for  q'  >  p'  there  is 
no  net  contribution  for  p'  >  k  because,  by  (3.6),  the  interaction  of  k', 
p* ,  q'  then  yields  conservative  transfer  of  energy  within  a  triad  of  three 
wavenumbers  all  >k.  Therefore  (3.7)  can  be  rewritten 

“  k  p '  +k ' 

n(k,t)  =  J  dk'J  dp'/  dq'T(k',p',q',t),  (3.8) 

k  0  (p£  =  larger  of  p’,  |p'-k'|). 

This  form  is  the  more  useful  one  for  computing  n(k,t).  Since  k'  is  always 

>k  and  p'  always  <k  in  (3,8),  contributions  from  very  small  or  very  large 
wavenumbers  always  involve  triads  with  a  big  ratio  of  largest  to  smallest 
wavenumber.  If  the  transfer  is  effectively  local  in  wavenumber,  the  inte¬ 
grations  in  (3.8)  converge  at  zero  and  infinity  in  such  fashion  that  only 
inertial-range  k',  p',  q'  contribute  if  k  lies  within  a  long  inertial  range. 
Equation  (3.7),  on  the  other  hand,  gives  a  nonlocal  representation  of  n(k,t)f 
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and  a  correct  result  is  obtained  from  it  for  inertial-range  k  only  if 
the  integrations  are  extended  to  include  the  entire  dissipation  range. 

Insight  about  the  structure  of  the  nonlinear  interaction  is  provided 
by  considering  in  detail  the  contributions  to  T(k,t)  from  interactions  in¬ 
volving  wavenumbers  very  much  larger  or  very  much  smaller  than  k.  Heis- 
,  3 

enberg  proposed  that  the  effect  of  higher  wavenumbers  on  T(k,t)  could 
be  represented  by  an  eddy  viscosity.  The  present  closure  approximation 
supports  this  hypothesis  for  wavenumbers  much  larger  than  k  and  yields  an 
analytical  expression  for  the  eddy  viscosity  (different  in  form  from 
Heisenberg's). 

A  Lagrangian  time-correlation  function  for  mode  k  can  be  defined  by 

R(k;t|r)  =  U(k;t|r)/[U(k;t|t)U(k;r|r)]1/2.  (3.9) 

To  evaluate  the  asymptotic  eddy  viscosity,  assume  first  that  both  E(k,t) 
and  the  characteristic  decay  times  t-r  of  G(k;t|r)  and  R(k;t|r)  decrease 
as  k  rises.  If  these  conditions  are  satisfied  with  sufficient  strength 
(to  be  checked  later),  the  first  term  on  the  right-hand  side  of  (2.6)  gives 
a  negligible  part  of  the  contribution  to  T(k,t)  from  wavenumbers  »k.  The 
reason  is  that  p  »  k  implies  q  >>  k,  and  then  U(k;t |s)U(q;t |s)  » 
U(p;t|s)U(q;t|s).  Furthermore,  U(k;t|s)»  U(k;t|t)  in  the  second  term, 
for  s  which  make  appreciable  contributions.  There  results 

®  q+k  t 

T>q'^k»t^  =  “  2E(k,t)/  ^ dqj  dpB  /  G(p;t | s)U(q;t | s)ds  (3.10) 

q'  q-k  t0 

for  the  contribution  to  T(k,t)  from  all  q  >  q'  >>  k. 
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The  integration  in  (3.10)  is  symmetric  in  p  and  q,  to  leading  order 
in  k/q.  To  leading  order  in  k/q, 

i(B.  -  B.  )  *  -  irq3(u  -  u3),  (3.11) 

2  kpq  kqp 

where  u  =  (p-q)A.  The  trigonometric  identity 

2AkPq  5  Bkpq  +  Bkqp  =  wkM(1  ’  ***  "  (3*12) 

yields  that  to  leading  order 

1  (B,  +  B,  )  »  ^  Aq2(l  +  u2  -  2u4).  (3.13) 

2  kpq  kqp  2  ^ 

If  the  G  and  U  factors  in  (3.10)  are  expanded  in  Taylor  series  about  u  =  0, 


these  expressions  lead  to 


where 


and 


T  (k,t)  =  -  2v  (t)k2E(k,t)  (q  »  k),  (3.14) 

>q  q 

00  “t  r\ 

v  (t)  =  —■  /  [7/  G(p;t|s)U(p;t|s)ds  -  pA(p,t)]p  dp  (3.15) 

<1  to 

.  t 

A(q,t)  =  {3^  J  [G(pjt|s)U(q;t|s)  -  G(q;t | s)U(p;t | s)]ds}p=q.  (3.16) 
t0 


Note  that  the  eddy  viscosity  v^(t)  is  independent  of  k  and  therefore  is  a 
true  counterpart  of  molecular  viscosity. 

Now  define  T<(^(k,t)  as  the  total  contribution  to  T(k,t)  due  to  all 
interactions  involving  a  wavenumber  less  than  q.  According  to  (3.14),  the 
loss  from  wavenumber  k  due  to  interaction  with  much  higher  wavenumbers  is 
proportional  to  the  mean-square  shear  or  vorticity  associated  with  k.  The 
interactions  are  conservative.  Therefore  T<(|k,t)  should  contain  a  contri¬ 
bution  that  involves  v^(t)  and  the  mean-square  vorticity  in  wavenumbers 
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less  than  q,  if  q  <<  k.  Actually,  the  complete  expression  for  T<(^(k,t) 
cannot  be  deduced  from  v^(t)  because  the  interaction  is  among  triads 
rather  than  pairs  of  wavenumbers.  The  eddy  viscosity  represents  the  net 
loss  of  energy  from  low  to  pairs  of  high  wavenumbers.  In  addition,  the 
interaction  with  a  low  wavenumber  induces  a  transfer  of  energy  between 
the  two  high  wavenumbers,  and  these  two  effects  turn  out  to  be  of  com¬ 
parable  importance  in  the  overall  energy  transfer  at  the  high  wavenumbers. 


The  complete  asymptotic  expression  for  T<(^(k,t)  is 
,  q*  ,  k+q  t 

T<Q,(k,t)  =  2k  /  E(q,t)q"  dqj  dp/  [2A,  G(k;t|s)U(p;t|s) 
q  0  k-q  tn  *P<1 


-  BkpqG(p;t|s)U(k;t|s)]ds  (k»q'),  (3.17) 

which  is  obtained  from  (2.6)  by  neglecting  the  term  involving  U(k;t|s) 
for  p  <  q'  and  assuming  that  the  characteristic  times  of  mode  q  are  long 
compared  to  those  of  k  and  p.  The  factor  2  in  the  first  term  on  the  right- 
hand  side  of  (3.17)  comes  from  including  the  contributions  of  both  q  <  q' 
and  p  <  q'  in  the  first  term  on  the  right-hand  side  of  (2.6).  There  is  a 
cancellation  between  the  A  and  B  terms  in  (3,17)  to  leading  order  in  p-k. 

The  remainder  can  be  evaluated  conveniently  by  considering  symmetrical 
and  antisymmetrical  parts  as  in  obtaining  (3.15).  The  final  result  can  be 
written 


(k  >>  q). 


(3.18) 
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The  three  terms  in  (3.18)  are  most  easily  interpreted  by  integrating 
over  k  to  form  the  corresponding  contribution  to  n(k,t): 

II  (k,t)  =  2f  E(p,t)p2dp[v.  (t)  +  ts-  k4A(k,t) 

<q  q  K  XD 

-21  k3/*  G(k;t|s)U(k;t|s)ds]  (k  »  q).  (3.19) 

*0 


The  first  term  on  the  right  clearly  corresponds  to  the  loss  from  low 
wavenumbers  given  by  (3.14).  If  the  modes  ■*'k  were  in  equipartition 
equilibrium  among  themselves,  they  would  satisfy  the  relations1 

G(k;t|s)  =  R(k;t|s),  3U(k;t | t)/Sk  =  0,  (3.20) 


whence  A(k,t)  would  vanish.  On  the  other  hand,  if  U(k;t|t)  is  a  rapidly 
decreasing  function  of  k,  as  in  the  inertial  range,  A(k,t)  is  positive. 
Its  appearance  in  (3.19)  then  represents  a  net  transfer  of  energy  from 
modes  in  the  range  (k-q,  k)  to  modes  >k  due  to  straining  by  modes  <q 
but  without  net  loss  from  the  modes  <q.  The  last  term  in  (3.19)  survives 
even  if  the  modes  ~k  are  in  equipartition.  Equation  (2.6)  is  easily 
shown  to  give  zero  transfer  among  a  triad  all  three  of  whose  members  are 
in  equipartition  equilibrium.  This  may  be  called  the  detailed  balance 
property.  The  presence  of  the  last  term  in  (3.19)  illustrates  that 
detailed  balance  between  only  two  members  of  a  conservative  interacting 
triad  of  modes  is  an  invalid  concept. 
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When  viscosity  effects  are  strong,  the  characteristic  fall-off  time  of 

G(q;t|s)  is  (vq2)”1,  and  the  result  is  a  value  of  which  is  smaller 

than  Heisenberg's  and  inversely  proportional  to  the  molecular  viscosity. 

The  absence  of  an  interference  between  molecular  and  eddy  viscosities  in 

Heisenberg's  formula  leads  to  the  prediction  of  a  power-law  dissipation- 

range  spectrum,  which  implies  non-existence  of  high-order  spatial  deriva- 

4  5 

tives  of  the  velocity  field.  ' 

The  asymptotic  formula  for  T<(^(k,t)  indicates  that  energy  transfer 
associated  with  widely  separated  wavenumbers  is  proportional  to  the  mean- 
square  vorticity  in  the  lower  wavenumbers.  It  is  important  in  understand¬ 
ing  the  structure  of  (2.6)  to  remember  that  this  result  involves  a  near 
cancellation  between  the  two  terms  on  the  right-hand  side.  These  terms 
separately  give  contributions  proportional  to  the  energy  in  the  wavenumbers 
<q,  not  to  the  mean-square  vorticity.  The  term  containing  G(k;t|s)  is  an 
input  term  giving  a  positive  contribution  to  T(k,t),  while  the  term  con¬ 
taining  U(k;t|s)  is  an  output  term  representing  a  drain  of  energy  from 
mode  k  to  other  modes.  The  input  and  output  contributions  from  low  wave- 
numbers  are  proportional  to  energy  because  they  represent  convection  as 
well  as  straining.  Convection  of  high-wavenumber  structures  by  strongly 
excited  low-wavenumber  velocity  components  implies;  a  rapid  change  of  phase 
of  the  high- wavenumber  Fourier  amplitudes;  that  is  to  say,  a  rapid  exchange 
of  energy  between  sine  and  cosine  components  of  the  high  wavenumbers.  This 
exchange  is  represented  in  (2.6)  by  the  large,  cancelling  input  and  output 
contributions.  The  net  contribution  T<(^(k,t)  represents  the  effect  of 
straining  alone. 


-16- 


Since  the  output  integral  in  (2.6)  involves  U(k;t|s),  it  gives  a 
negative  -contribution  to  T(k,t)  roughly  proportional  to  E(k,t).  As  a 
result,  the  general  behavior  of  (2.6)  is  to  give  energy  flow  out  of 
strongly  excited  modes  into  weakly  excited  modes.  In  equipartition, 
the  input  and  output  terms  cancel  exactly.  This  can  be  seen  immediately 
by  using  (3.20)  in  (2.6),  as  noted  earlier.  Out  of  equipartition,  the 
energy  flow  is  from  low  to  high  wavenumbers  only  if  U(k;t|t)  decreases 
as  k  increases. 


4.  THE  TIME-CORRELATION  AND  RESPONSE  EQUATIONS 

The  multiplicity  of  terms  in  (2.7)  and  (2.8)  arises  from  the  peculiar 
way  in  which  Lagrangian  dynamics  are  represented  through  (1.2).  A  more 
standard  Lagrangian  analysis®  is  based  on  the  equation  for  the  particle 
acceleration  3u(x,tQ |r)/3r.  This  treatment  leads  to  trouble  in  turbulence 
theory  because  the  expressions  for  viscous  and  pressure  forces  involve 
awkward  nonlinearity,  which  so  far  has  precluded  approximate  statistical 
equations  that  preserve  incompressibility  and  conserve  energy.  In  the 
present  approach,  this  trouble  is  sidestepped  by  avoiding  direct  calcula¬ 
tion  of  3u(x,t|r)/3r.  Instead,  the  basic  equations  are  for  3u(x,t |t)/3t 
and  3u(x,t |r)/3t.  The  r  dependence  is  determined  by  the  difference 
between  the  behaviors  along  the  appropriate  paths  of  integration  of  these 
two  equations.  The  different  approaches  are  illustrated  in  Fig.  1. 

Suppose  that  the  initial  Eulerian  field  u(x,tQ)  =  u(x,t0|tQ)  is  given  and 
u(x,t0|r)  is  desired.  The  straightforward  procedure  is  to  integrate 
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3u(x,tQ|r)/3r  along  the  path  AB.  In  the  present  treatment,  the  Navier- 
Stokes  equation  is  integrated  along  AC,  yielding  u(x,t|t)  for  tQ  <  t  <  r. 

Then  (1.2)  is  integrated  backward  in  t  along  CB.  Thus  u(x,tQ |r)-u(x,tQ |tQ) , 
the  total  change  in  the  particle  velocity  during  the  interval  (tQ,  r), 
appears  as  the  difference  between  the  two  forward  integrations  AC  and  BC. 

If  there  are  no  viscous,  pressure,  or  external  forces,  the  two  integrations 
are  identical,  so  that  u(x,tQ|r)-u(x,t0 | tQ)  =  0. 

Since  (x,t)  is  laboratory  spacetime,  this  amounts  to  describing 
Lagrangian  behavior  by  equations  which  are  Eulerian  in  character,  a  fact 
pointed  out  to  the  author  by  Prof.  S.  Corrsin.  The  advantage  is  that 
conservation  and  invariance  properties  appear  directly  and  are  easily 
incorporated  in  the  final  closure  equations.  The  penalty  is  that  the 
Lagrangian  statistical  functions  are  determined  as  algebraic  sums  of  a 
multiplicity  of  intricately  counterbalancing  terms. 

The  meaning  of  the  various  terms  in  (2.7)  can  be  uncovered  by  tracing 
how  they  arise  in  the  perturbation  analysis"'"  through  which  the  equation 
was  derived.  The  v  term  has  an  obvious  significance:  simple  viscous  decay 
of  Fourier  amplitude  of  the  Eulerian  field  at  time  t.  The  first  term  on 
the  right-hand  side  can  be  interpreted  by  considering  a  hypothetical  situation. 
Imagine  that  a  velocity  field  with  an  isotropic,  multivariate-Gaussian 
ensemble  distribution  is  set  up  at  time  r  and  then  forces  are  somehow 
applied  so  that  the  Eulerian  field  is  frozen  thereafter,  exhibiting  no 
time -dependence  at  all.  The  fluid  particles  move  along  the  streamlines  of 
the  spatially  varying  Eulerian  field,  and  the  Lagrangian  velocities  change 
with  time.  If  the  derivation  of  (2.7)  is  repeated  for  this  situation. 


the  only  surviving  term  is  the  first  term  on  the  right-hand  side.  This 
term  represents  a  decrement  in  the  correlation  U(k;t|r)  due  to  the  labeling 
transformation  (transformation  from  Eulerian  to  Lagrangian  coordinates): 

The  particle  at  x  at  time  t  was  elsewhere  at  time  r,  so  that  u(x,t|r)  is 
really  the  Eulerian  velocity  at  a  point  other  than  x  and  is  partly  decorrelated 
with  u(x,t|t). 

In  the  actual  case,  where  the  Eulerian  velocity  is  not  static,  this 
correlation  loss  is  partly  compensated  by  a  correlation  of  the  relabeling 
with  the  evolution  of  the  Eulerian  field  during  the  same  period:  The 
substantial  derivative  [3/3t  +  u(x,t)*2]  which  determines  the  relabeling 
also  appears  in  the  Navier-Stokes  equation.  The  second  term  on  the  right- 
hand  side  of  (2.7)  expresses  the  compensation.  The  first  and  second  terms 
together  express  net  effects  of  relabeling  during  the  period  (r,t).  Call 
them  the  Class  I  terms. 

The  remaining  terms  on  the  right-hand  side  of  (2.7)  all  express 

correlations  arising  from  the  distortion  of  the  Eulerian  field  over  its 

entire  history  from  tQ.  They  consist  of  Class  II,  with  upp  r  limit  r,  and 

Class  III,  with  upper  limit  t.  Each  of  these  classes  is  itself  a  sum  of 

compensating  terms.  The  algebraic  sum  of  the  terms  within  each  class 

provides  a  cancellation  which  eliminates  convection  effects  in  the  Eulerian 

evolution;  only  the  difference  between  the  histories  at  different  space 

points  makes  up  the  net  contribution.  The  integrals  of  Class  II  and  Class 
3 

III  over  d  k  separately  vanish  for  all  t  and  r.  This  expresses  that  past 
nonlinear  processes  in  the  Eulerian  field  have  no  direct  effect  on 
3(u(x,t|t)*u(x,t|r)^/3t.  The  latter  derivative  is  determined  wholly  by 
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the  labeling  transformation  in  the  interval  (r,t).  The  vanishing  of  the 

3 

integrals  over  d  k  can  be  verified  by  using  (2.9),  (3.4),  and  some  changes 
of  variable. 

It  may  be  asked  why  the  net  effect  of  past  distortion  of  the  Eulerian 
field  does  not  appear  simply  as  an  integral  over  the  interval  (r,t)  instead 
of  as  the  algebraic  sum  of  Class  II  integrals  over  (tQ,r)  and  Class  III 
integrals  over  (tQ,t).  The  reason  is  that  distortion  during  (r,t)  is 
affected  by  higher  correlations  built  up  over  the  entire  history  of  the 
dynamical  interaction  as  well  as  by  the  value  of  the  velocity  covariance 
at  time  r.  It  will  be  seen  below  that,  in  general,  the  Class  II  and 
Class  III  terms  do  not  cancel  each  other  if  t  =  r.  They  express  effects 
that  exist  even  if  no  labeling  transformation  has  taken  place. 

A  similar  analysis  may  be  made  for  the  response  equation  (2.8).  The 
first  two  terms  on  the  right-hand  side  arise  from  relabeling  associated 
with  the  unperturbed  velocity  field;  that  is,  relabeling  of  the  perturba¬ 
tion  according  to  the  transformation  that  would  take  place  in  absence  of 
perturbation.  They  are  analogous  to  the  Class  I  terms  in  (2.7).  In 
addition,  the  applied  perturbation  at  time  r,  to  which  G(k;t|r)  gives  the 

average  response,  induces  a  perturbation  in  the  relabeling  transformation. 

term 

This  is  expressed  by  the  [D^  -  Bkpq^on  the  right-hand  side  of  (2.8). 

The  remaining  two  terms  represent  effects  from  the  perturbation  in  the 
nonlinear  evolution  of  the  Eulerian  field. 

Hopefully,  the  interpretations  above  help  a  little  to  make  (2.7)  and 
(2.8)  intelligible.  However,  no  clean-cut  separation  is  possible  among 
the  various  effects  cited.  All  the  terms  together  determine  the  values  of 
the  G  and  U  functions  which  appear  in  any  one  term. 
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It  was  stated  in  Sec.  1  that  ^u(x,t |t) *u(x' ,t |r)^  was  not  symmetric 
about  r  =  t  in  stationary ,  homogeneous  turbulence ,  except  if  x  =  x' .  This 
can  be  investigated  with  (2.7).  First,  an  isotropic  energy  source  is  needed 
in  order  that  stationarity  and  isotropy  be  mutually  consistent.  The  easiest 
way  is  to  add  a  negative  damping  term  to  the  Navier-Stokes  equation  so 
that,  after  Fourier  transformation, 

vk2  +  vk2  -  u(k),  (4.1) 

2 

where  u(k)  is  positive,  and  >vk  for  some  k.  Let  u(k)  vanish  outside  the 

energy-containing  range.  This  input  mechanism  may  be  regarded  as  an 

idealization  of  the  mechanism  by  which  real  turbulence  draws  energy  from 

an  overall  shear.  The  effect  of  the  input  on  (2.6)-(2.8)  is  simply  to 
2  2 

replace  vk  by  vk  -  u(k)  everywhere.  To  describe  the  stationary  state, 
the  limit  tQ  is  taken. 

Now  consider  (2.7)  for  r  =  t.  The  two  D'  terms  cancel  exactly,  and, 
by  (2.9)  and  the  symmetry  of  the  integrand  factors  in  p  and  q,  the  coeffi¬ 
cient  i°  the  sixth  term  can  be  replaced  by  Bj^.  The  first  two  terms 

on  the  right-hand  side  vanish.  The  result  is  that  the  right-hand  side  is 
identical  with  the  right-hand  side  of  (2.6).  Now  subtract  (2.6)  from  twice 
(2.7)  and  multiply  by  2wk2.  In  the  stationary  state,  3U(k;t|t)/3t  =  0. 
Therefore,  by  the  definitions  of  T(k,t)  and  R(k;t|r),  the  result  is 

E(k,t)[3R(kjt|r)/3t]r_t  =  yT(k,t), 
or,  by  stationarity, 

E(k,t)[3R(k;t|r)/3r]r_t  =  -  ~  T(k,t). 


(4.2) 


Fourier  transformation  of  (4.2)  yields 


[3Ui.(x,t|t;x+jj,t|r)/3r]r=t  =  -/  T(k,t)S1j^-ty-'  dk.  (4.3) 

For  y  =  0,  (4.3)  gives  zero  slope  at  r  =  t,  by  (3.3).  For  y  intermediate 
between  energy-containing  and  dissipation  scales,  the  sin(ky)/(ky)  factor 
depresses  the  relative  contribution  of  wavenumbers  receiving  energy,  so 

that  the  right-hand  side  of  (4.3)  is  positive.  This  agrees  with  the 

.  .  ,  .  ,  „  wavenumbers 

discussion  in  Sec.  1.  Equation  (4.3)  shows  thatAwhere  energy  is  fed  in 

[negative  T(k,t)]  contribute  acceleration  positively  correlated  with 
particle  velocity,  while  wavenumbers  dissipating  energy  contribute  nega¬ 
tively  correlated  acceleration. 

Equation  (4.2)  actually  is  an  exact  result  of  stationarity,  independ¬ 
ent  of  the  closure  approximation  and  of  the  nature  of  the  input  forces. 

If  an  isotropic,  statistically  stationary  forcing  term  of  any  kind  is 
added  to  the  Navier-Stokes  equation,  then  (4.2)  follows  immediately  from 
the  identity  of  the  triple  moments  and  input  terms  in  the  exact  equations 
for  3U(k;t|t)/3t  and  [3U(k;t|r)/3t]r_t  that  can  be  formed  from  (1.2)  and 
the  Navier-Stokes  equation. 

Evenness  of  U^. (x,t | t;x,t |r)  as  a  function  of  r-t  requires  that  all 
the  higher  odd  r  derivatives  vanish  at  r  =  t,  as  well  as  the  first.  The 
behavior  of  the  higher  derivatives  has  not  been  determined  for  either 
(2.6)-(2.8)  or  the  unabridged  LHDI  equations.  This  question  is  part  of 
the  general  problem  of  integrability  and  realizability  of  the  LHDI  equations 
discussed  in  Ref.  1. 

Further  insight  about  the  time-correlation  and  response  equations  is 
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provided  by  investigating  the  contributions  from  interactions  with  much 
larger  and  much  smaller  wavenumbers.  It  is  interesting  first  to  examine 
the  contribution  from  very  local  interactions,  p  ss  k,  q  as  k.  Suppose  that 
all  triad  interactions  are  absent  from  (2,6)-(2.8)  except  those  for  which 

|p  -  k|  <  6,  |q  -  k |  <6,  6  <<  k.  (4.4) 

Let  6  be  small  enough  that  p  and  q  may  be  replaced  by  k  everywhere  in  the 
integrands .  Then 

B,  as  D.  ®  Di  as  I-  irk3,  C.  *  £  irk3.  (4.5) 

kpq  kpq  kpq  8  *  kpq  4 

The  right-hand  side  of  (2.6)  vanishes  to  this  approximation;  because  of 
symmetry,  very  local  interactions  produce  no  net  energy  transfer.  All 
terms  but  the  first  two  in  the  right-hand  sides  of  (2.7)  and  (2.8)  either 
cancel  in  pairs  or  vanish.  When  the  integrations  over  p  and  q  are  per¬ 
formed,  the  net  contributions  from  the  first  two  terms  yield 

(3/3t  +  vk2)U(k;t|r)  =  -  k362U(k;t  |r)J  U(k;t|s)ds,  (4.6) 


(3/3t  +  vk2)G(k;t|r)  =  -  k3fi2G(k;t  |r)J  U(k;t|s)ds, 

r 

G(k;r |r)  =  1.  ( 4 

The  exact  solution  to  (4.6)  and  (4.7)  for  v  =  0  can  be  obtained  by 
the  trial  substitution 

U(k;t|r)  =  U(k)G(k;t|r)  =  g(t-r)U(k),  ( 4 

which  reduces  (4.6)  and  (4.7)  to  the  single  equation 


3g(t)/3t  =  -  a  g(t)J  g(s)ds,  g(0)  =  1, 
0 
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where  a2  =  ( 3ir/2)k^62U(k) .  The  solution  of  (4.9)  is 

g(t)  =  4e‘att/2/(l  +  e‘at/2)2,  (4.10) 

which  is  an  even  function  of  t  and  has  the  asymptotic  behavior 

g(t)  »  1  -  j  (at)2  (at<<l),  g(t)  »  4e“atl,/2  (at>>l).  (4.11) 

This  solution  suggests  that  interactions  local  in  wavenumber  tend  to  pro¬ 
duce  an  eventual  rapid  decay  of  U(k;t|r)  and  G(k;t|r)  at  large  t-r.  Setting 

6  =  k  suggests  that  the  characteristic  decay  frequency  due  to  interactions 

2 

with  wavenumbers  of  order  k  is  crudely  kv^,  where  v^  =  kE(k,t). 

To  obtain  the  contribution  to  (2.7)  from  all  interactions  involving 
wavenumbers  >q  >>  k,  make  the  same  assumptions  used  in  obtaining  (3.10). 

The  asymptotic  contribution  is 

«  q'+k  t 

[3U(k;t |r)/3t]  =  -U(k;t | r) /  dq'J  dp|-k  , /  U(q’;t|s)ds 

**  q  q'-k  ”  r 

t 

+  B,  ,f  G(p;t)s)U(q' ;t|s)ds]  (q>>k).  (4.12) 

KP<1  4. 

0 

All  the  other  terms  on  the  right-hand  side  are  suppressed  either  because 
they  involve  U(p)U(q)  or  because  they  involve  at  least  one  function  for 
p  or  q  with  a  time  separation  of  order  t-r;  (4,12)  assumes  that  t-r  is 
large  compared  to  the  characteristic  decay  times  of  the  modes  >q.  In  the 
second  term  of  (4.12),  U(k;r|s)  «  U(k;r|t)  =  U(k;t|r)  has  been  used,  because 
only  s  *  t  can  contribute.  The  wavenumber  integrations  in  this  term  were 
performed  in  evaluating  (3.10).  The  integration  in  the  first  term  is  easily 
performed,  with  the  result 
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[3U(k;t|r)/3t]>q  =  (t)U(k;t  |r)  (q  >>  k). 


v'(t)  =  v  (t)  +  J  dq'J  ( q * ) 2U( q* ;t|s)ds. 


(4.13) 


In  the  expression  for  v^(t),  the  lower  limit  r  has  been  changed  to  tQ. 
This  is  permissible  because  all  appreciable  contribution  comes  from 
t-s  <<  t-r.  Similar  analysis  gives 


[3G(k;t|r)/3t]>(^  =  -  v^(t)G(k;t  |r)  (q  >>  k)  (4.14) 

for  the  contribution  to  the  right-hand  side  of  (2.8)  from  wavenumbers  aq, 
for  t-r  large  compared  to  the  characteristic  times  of  these  wavenumbers. 

Comparison  of  (4.13)  and  (4.14)  with  the  viscous  terms  in  (2.7)  and 
(2.8)  shews  that  the  high-wavenumber  interactions  act  like  an  eddy  viscosity, 
but  not  the  same  eddy  viscosity  which  appears  in  the  energy-transfer  con¬ 
tribution  (3.14).  The  difference  may  be  understood  as  follows.  The 
energy  transfer  among  different  wavenumbers  in  isotropic  turbulence  can  be 
ascribed  to  advective  straining.  The  pressure  term  makes  no  direct 
contribution.  The  pressure  term  also  makes  no  direct  contribution  to  the 
evolution  of  the  Eulerian  time  covariance  ^u(x,t)*u(x' ,t' )).  The  exact 
equation  for  this  covariance  is 

( 3/3t  -v7^)(u(x,t ) *u(x' ,t'  ))  =  -  (u(x,t) *Vxu(x,t) *u(x' ,t'£> 

-  ^u(x' ,t,),Vp(x,t)^>>  (4.15) 


and  the  term  containing  the  pressure  p(x,t)  vanishes  by  incompressibility 
and  homogeneity. 

Just  the  reverse  is  true  for  the  Lagrangian  time  covariance.  Advection 
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produces  no  direct  change  in  particle  velocity  and  therefore  makes  no 

direct  contribution  to  decay  of  the  Lagrangian  covariance.  The  pressure 

fluctuations,  on  the  other  hand,  produce  accelerations  which  decorrelate 

the  Lagrangian  velocity.  Crudely  speaking,  v^(t)  may  be  called  the  eddy 

viscosity  due  to  advective  straining  and  v^(t),  the  eddy  viscosity  due 

to  pressure  fluctuations.  The  separation  is  not  clean  because  the  advec- 

tion  and  pressure  fluctuations  react  on  each  other. 

The  difference  v^(t)  -  v^(t)  is  a  positive  integral.  This  provides  a 

crude  suggestion  that  the  Lagrangian  correlation  time  for  turbulence  may 

be  shorter  than  the  Eulerian  correlation  time.  Some  simple  physical 

7 

reasoning  leads  to  the  same  suggestion.  The  argument  is  that  even  in  a 
frozen  Eulerian  velocity  field  (infinite  Eulerian  correlation  time), 

Lagrangian  covariances  decay  with  difference  time  because  given  fluid 
elements  migrate  through  the  velocity  field.  This  question  will  be  dis¬ 
cussed  further  in  Sec.  6,  after  the  quantitative  inertial-range  results 
are  available. 

Now  consider  the  contributions  to  (2.7)  and  (2.8)  from  interactions 
with  wavenumbers  Sq’  <<  k.  For  t-r  the  order  of  the  correlation  and  response 
times  of  mode  k,  the  situation  is  like  that  for  energy  transfer.  The  contri¬ 
butions  to  the  various  terms  on  the  right-hand  sides  cancel  in  such  fashion 

that  what  remains  is  proportional  to  the  mean-square  vorticity  in  the  wave- 

terms 

numbers  <q' .  In  (2,7),  the  Class  I,  Class  II,  and  Class  IIIAcancel  separately. 
In  (2.8),  there  are  cancellations  between  the  first  two  terms,  within  the 
third  term,  and  between  the  last  two  terms. 

The  behavior  at  large  t-r  is  more  complicated.  The  analysis  above 
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suggests  that  if  k  interacted  only  with  wavenumbers  of  order  k  and  larger, 
the  functions  U(k;t|r)  and  G(k;t|r)  would  fall  off  exponentially  for  large 
t-r.  What  actually  happens  is  that  interactions  with  low  wavenumbers 
dominate  the  long-time  behavior  and  cause  these  functions  to  fall  off 
more  slowly.  Discussion  of  this  phenomenon  is  postponed  to  Sec.  6. 


5 .  INERTIAL-RANGE  EQUATIONS 

Kolmogorov's  inertial-range  law  is 

E(k)  =  Ce2/3k-5/3,  (5.1) 

where  e  is  the  rate  of  energy  dissipation  per  unit  mass  and  C  is  a  universal 
constant.  The  time-dependence  is  omitted  in  (5.1)  because  the  relevant 
dynamical  times  in  the  inertial  range  are  assumed  to  be  short  compared  to 
the  decay  time  of  the  total  turbulent  energy,  so  that  the  inertial  range  is 
nearly  statistically  stationary.  Kolmogorov's  hypotheses  also  yield 
inertial-range  forms  for  G(k;t|r)  and  R(k;t|r).  The  requirement  that  they 
depend  only  on  e,  k,  and  t-r  imposes  on  these  functions  the  forms 

G(k;t|r)  =  G(c1/3k2/3x) ,  R(k;t|r)  =  R( e1/3k2/3x) ,  x  =  t-r.  (5.2) 

Reference  1  argued  that  the  LHDI  equations  led  to  a  Kolmogorov  inertial 
range  because  they  gave  the  required  localness  of  the  dynamical  interaction 
in  wavenumber.  This  can  be  verified  for  the  abridged  LHDI  equations  by 
substituting  (5.1)  and  (5.2)  into  (2.6)-(2.8)  and  checking  that  the  inte¬ 
grations  over  p  and  q  converge  properly  at  zero  and  infinite  wavenumbers. 

The  equations  for  G(x)  and  R(x)  which  result  from  the  substitution  can  be 
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written  most  usefully  by  introducing  rescaled  functions 

G(t)  =  G( YT) ,  R(t)  =  R(yO,  y  =  (2*/C)1/2,  (5.3) 

1/3  2/3 

where,  for  convenience,  t  now  denotes  the  nondimensional  argument  e  k  x 
of  (5.2).  Then  C  drops  out  of  (2.7)  and  (2.8)  leaving 

T  T 

dG(x)/dx  =  JJ  dvdw{-C  G(x)/  R(w2/3s)ds  +  D  G(v2/3t)J  R(w2/3s)ds 

4  4  A  VW  4  VW  '  q 

-  [Duv  -  Buv]G(v2/3t)R(w2/3t)J  G(s)ds 
-  /  [B  G"(v2^3s)G"(t-s)  -  D'  G(s)G"(v2^3(T-s))]R(w2,/3s)ds} 

'  Q  VW  VW 

(t  >  0),  G(0)  =  1,  (5.4) 

T  T 

dR(x)/dx  =  //. dvdw{-C  R(t)/  R(w2//3s)ds  +  D  v  11//3R(v2/,3t)J  R(w2^3s)ds 

n  vw  o  w  0 

00 

+  /  [B  v_11^3G^(s-t)R(v2'/3s)  -  D'  G(v2/,3(s-x))R(s)]R(w2/,3s)ds 

4  VW  VW 

T 

00 

-  f  [B  G’(v2/,3s)R(  | s-t  | )  -  D  v"11^3g’(s)R’(v2'/3  |  s-t  |  )]R(w2^3s)ds 

4  Q  VW  VW 

oo 

-  /  [B  -  D’  ]G’(v2/'3s)R(s)R(w2/,3|s-t|  )ds}  (x  z  0),  R(0)  =  1.  (5.5) 

4  0  VW  VW 

Here  v  =  p/k,  w  =  q/k,  and  JJ  now  denotes  integration  over  all  of  the 
(v,w)  plane  such  that  (l,v,w)  can  form  the  sides  of  a  triangle.  In  accord 
with  the  assumption  of  an  asymptotically  stationary  inertial  range,  the 
limit  tQ  -*■  -»  has  been  taken.  In  addition,  some  simple  transformations  of 
time  variables  have  been  made  and  the  symmetry  property  (2.5)  has  been  used. 
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The  coefficients  in  (5.4),  (5.5)  can  be  written 

n  2  -8/3,  3.  „  -8/3,,  2, 

Bvw  =  v  w  (xy  +  z  ),  Cw  =  7TVW  (1  -  y  ), 

Dvw  =  7  tivw"8//3(1  “  xyz  -  2y2z2  +  z2  -  y2), 

Dvw  =  J  1Tv3w”8^3(l  -  xyz  -  2x2z2  +  z2  -  x2),  (5.6) 

where  x,  y,  z  are  now  the  cosines  of  the  interior  angles  opposite  the 
triangle  sides  1,  v,  w  respectively. 

Substitution  of  (5.1)  and  (5.2)  into  (2.6)  makes  T(k,t)  vanish 
identically  by  virtue  of  the  conservation  identity  (3.4).  The  energy 
gained  by  k  from  lower  modes  is  exactly  balanced  by  a  drain  to  higher  modes. 
This  formal  result  is  a  necessary  consistency  requirement  for  an  inertial 
range,  but  it  says  nothing  about  the  localness  of  the  energy  transfer. 
Transfer  between  widely  separated  wavenumbers  can  be  investigated  by  sub¬ 
stituting  (5.1)  and  (5.2)  into  (3.15)  and  (3.19).  The  result  is  that  both 

the  loss  from  modes  <q  due  to  interactions  with  modes  2k  and  the  gain  in 

4/3 

modes  >k  due  to  interactions  with  modes  <q  are  “(q/k)  for  k  >>  q,  pro¬ 
vided  that  G(t)  and  R(t)  are  properly  integrable  over  t.  Discussion  of  the 
behavior  of  the  latter  functions,  and  the  check  of  localness  for  (5.4)  and 
(5.5),  are  discussed  later  in  this  Section  and  in  Sec.  6. 

The  Kolmogorov  constant  C  can  be  determined  by  evaluating  (3.8).  In 
the  inertial  range,  JI(k,t)  must  be  independent  of  k  and  equal  to  e,  since 
there  is  negligible  dissipation  and  the  energy  cascade  is  conservative. 

This  implies  the  similarity  law 

-3 

T(ak,  ap,  aq,  t)  =  a  T(k,  p,  q,  t), 


(5.7) 
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where  a  is  such  that  all  wavenumbers  involved  are  in  the  inertial  range. 
Equation  (5.7)  can  be  verified  for  the  abridged  LHDI  transfer  function  by 
inserting  (5.1)  and  (5.2).  It  leads  to  an  important  simplification  of  (3.8). 
Write  b  =  p'k/k',  q  =  q'k/k'  so  that  the  expression  for  e  =  n(k,t)  becomes 
»  k2/k'  b+k 

e  =/  (k/k'Jdk'J  db/  dqT(k,b,q)  ( bis  =  larger  of  b,  k-b),  (5.8) 

k  0  b#. 

by  (5.7).  The  time  argument  in  (5.8)  and  later  equations  is  suppressed  as 
in  (5.1).  Now  set  p  =  k2/k'  in  (5.8)  and  find 
k  p  b+k 

e  =  /  (k/p)dpf  db/  dqT(k,b,q).  (5.9) 

0  0  b* 

g 

Finally,  integrate  by  parts  and  obtain 
k  p+k 

e  =  kf  dp  ln(k/p)/  T(k,p,q)dq.  (5.10) 


If  (5.1)  and  (5.2)  are  substituted  into  (5.10),  some  simple  changes  of 
variable  yield  the  formula 


1+v 


(C/2ir)-^2  =  4 it /  dv  ln(l/v)J  dw/  ds[A  G(s)R(v2//3s)R(w2/,3s) 


v*  o 


-B  G(v2/3s)R(s)R(w2/3s)  -  B  G(w2/3s)R(s)R(v2/3s)],  (5.11) 

vw  wv 


where 


A  =  (B  +  B.  )(vw)”8//3, 
vw  lvw  lwv 


j  1 1  1 

vft  =  lv  -  2*  +  2 


It  should  be  remembered  that  although  T(k,t)  vanishes  identically  in 
the  inertial  range,  II(k,t)  defined  by  (3.7)  does  not  vanish,  because  the 
integrations  include  the  dissipation  range.  The  passage  from  (3.7)  to  (3.8) 
uses  the  conservation  properties  to  obtain  a  formula  for  C  which  involves 
only  inertial-range  functions. 
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Equation  (5.10)  can  be  cast  into  the  form 


e 


(5.12) 


where  H(o)(da/a)  is  defined  as  the  contribution  to  e  from  all  triad  inter¬ 
actions  such  that  the  ratio  of  maximum  to  minimum  wavenumber  lies  between  a 
and  a  t  da.  To  obtain  W(a),  write  (5.10)  as 

k  k  k  ktp 

e  =  kj  dp  ln(k/p)J  T(k,p,q)dq  +  kj  dp  ln(k/p )/  T(k,p,q)dq, 

0  Pa  Ok 


In  the  first  term  the  ratio  of  maximum  to  minimum  wavenumber  is  a  =  k/p 

and  in  the  second  term  it  is  a  =  q/p.  In  the  first  term  write  q  =  k/B 

and  in  the  second  write  p  =  k/B,  q  =  ak/B.  Then  some  manipulations  of 

the  limits  yields  e  in  the  form  (5.12)  with 
ft 

q  —1  ^  ®  q 

W(a)  =  k  [a-  lnaj  T(k,k/a,k/8)8  dB  +  a/  T(k, k/B, ak/B) B  InB  dfi], 

1  a-1 


*  .  rll  1  ,1,-1 

a  '  C'a  "  2  +  2  ]  * 


(5.13) 


which  is  independent  of  k,  by  (5.7).  Substitution  of  (5.1)  and  (5.2)  into 
(5.13)  yields 

Ha)  ■  a  1+1/3  a  >>  1,  (5.14) 

in  agreement  with  the  results  obtained  from  (3.15)  and  (3,19). 

For  values  of  t  small  enough  that  fT(t)  and  G(t)  are  still  0(1),  the 

convergence  properties  of  the  wavenumber  integrals  in  (5.4)  and  (5.5)  are 

like  those  of  the  energy  transfer.  Equations  (4.13)  and  (4.14)  yield 
-4/3 

contributions  «a  from  interactions  with  v  or  w  >  a,  a  >>  1.  Investi¬ 
gation  shows  that  the  contribution  from  interactions  involving  v  or  w  <a-1 
-4/3 

is  <*a  also.  The  behavior  at  large  t  will  be  discussed  in  Sec.  6. 
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6.  BASIC  NUMERICAL  RESULTS  AND  THEIR  INTERPRETATION 

Before  inertial-range  computations  were  attempted,  an  overall  consis¬ 
tency  check  of  the  closure  equations  was  performed  by  repeating  some  numerical 

studies  originally  carried  out  for  the  unaltered  Eulerian  direct-interaction 
g 

equations.  Inviscid  equipartition  solutions  were  investigated  first.  The 
viscosity  was  set  equal  to  zero  in  (2.6)-(2.8)  and  the  equations  truncated 
to  a  finite  wavenumber  range,  as  in  Ref.  9.  Initial  spectra  U(k;tQ|t0)  of 
several  forms  were  taken,  including  sharply  peaked  spectra  confined  to  a 
narrow  wavenumber  band.  The  system  was  found  to  approach  equipartition 
and  fluctuation-dissipation  equality  [Eq.(3.20)].  E(k,t)  stayed  positive 
for  all  k  and  t,  and  G(k;t|r)  and  R(k;t|r)  approached  zero  for  large  t-r. 

Some  of  the  finite-viscosity  decay  calculations  at  moderate  Reynolds 
number  described  in  Ref.  9  were  also  repeated.  The  spectral  predictions 
of  the  present  abridged  LHDI  equations  turned  out  to  be  qualitatively 
similar  to  those  of  the  unaltered  Eulerian  direct-interaction  equations 
at  the  Reynolds  numbers  investigated  (R^  &  40). 

The  unaltered  direct-interaction  equations  are  the  exact  description 
of  a  model  dynamical  system  and  are  guaranteed  self-consistent.  There  is 
no  such  guarantee  for  the  LHDI  equations.^"  The  numerical  results  just 
described,  and  those  to  follow,  suggest  that  the  LHDI  equations  actually 
are  a  more  stable  system  than  the  unaltered  direct-interaction  equations 
and  stay  farther  from  unphysical  behavior.  Functions  which  decay  with 
damped  oscillations  in  the  unaltered  approximation  instead  fall  smoothly  to 
zero  in  the  LHDI  and  abridged  LHDI  approximations.  A  possible  explanation 
lies  in  the  general  nature  of  the  alteration  prescription.  The  change  from 
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direct-interaction  to  LHDI  triple-moment  formulas  consists  of  replacing 

certain  labeling  times  s  <  t  by  t  in  the  integrals  over  the  history  of  the 

fluid.  This  appears  to  be  stabilizing  in  the  same  way  as  a  change  from 

forward  to  backward  differences  can  be  stabilizing  in  the  numerical 

integration  of  differential  equations. 

The  basic  inertial-range  equations  (5.4)  and  (5.5)  were  integrated 

by  a  modification  of  numerical  methods  developed  for  the  Eulerian  direct - 

g 

interaction  equations.  The  Kolmogorov  constant  was  then  computed  by 
evaluating  the  multiple  integral  (5.11)  in  such  a  way  that  the  calculation 
also  yielded  W(a)  without  extra  work.  The  results  for  G,  R,  W,  and 

Ka)  5  /  W(a')^- 

are  plotted  in  Figs.  2  and  3.  The  value  computed  for  the  Kolmogorov 

,  „  .  10 
constant  is 

C  =  1.77.  (6.1) 

A  description  of  the  numerical  procedures  and  estimates  of  computational 
errors  are  given  in  the  Appendix. 

The  following  values  for  integrals  over  R  and  G  were  obtained: 

00  00  oo 

I  =  /  G(s)ds  =  1.04,  I  =  /  R(s)ds  =  1.87,  I  =  /  G(s)R(s)ds  =  0.76, 
G  0  0  0 

CO 

I  =  /  [R(s)dG(s) /ds  -  G(s)dR(s) /ds]sds  =  -0.19.  (6.2) 

2  0 

These  numbers  yield  the  following  inertial-range  values  for  the  quantities 
which  appear  in  (3.15),  (3.18),  and  (4,13): 
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v.  =  i  [I  -  i  I  ]C  e1/,3k_4/3  =  0.235  e1/3k-4/3, 
k  b  1  5  2  ’ 


oo  t 

/  p2dp/  G(p;t|s)U(p;tjs)ds  =  —  CI1  e1/3k' 


=0.47  e  k 


J  p3A(p)dp  =  ~  [UI1  +  2I2]C  e1/3k“4/3  =  0.235  e1/3k“4/3, 
k 

|ik%(k)  =  [11^  +  2I2]C  e1/3k"4/3  =  0.314  e1/3k“4/3, 


f£k3J  G(k;t|s)U(k;t|s)ds  =  i_  e1/3k"4/3  =  0.090  E1/3k“4/3, 


v-  =vk+iCIR  e1/3k-4/3  =  1.89  c1/3k-4/3. 


(6.3) 


Figure  3  indicates  that  W(u)  peaks  at  a  =  2.  However,  it  also  indi¬ 
cates  that  only  about  20%  of  the  total  transfer  comes  from  a  <  2;  that  is, 

from  triad  interactions  in  which  the  largest  and  smallest  wavenumbers  differ 
by  a  factor  of  two  or  less.  The  tail  of  W(a)  carries  most  of  the  energy 
transfer,  so  that  the  local  cascade  appealed  to  by  Kolmogorov  seems  actually 
to  be  rather  diffuse.  This  may  mean  that  inertial  ranges  must  be  quite 

extensive  before  they  exhibit  nearly  asymptotic  dynamical  properties.  The 

-4/3 

numerical  results  support  the  a  dependence  of  W(a)  at  large  a  found 
analytically. 

The  numerical  results  (6.3)  for  effective  eddy  viscosities  and  other 
functions  relevant  to  distant  interactions  exHbit  some  interesting  features. 
Substitution  of  these  results  into  (3.19)  indicates  that  the  eddy-viscosity 
loss  from  low  wavenumbers  (the  first  term  on  the  right-hand  side)  represents 
slightly  more  than  50%  of  the  energy  gain  of  the  high  wavenumbers.  An 
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almost  equal  gain  in  the  wavenumbers  >k  comes  from  the  straining  action 
of  the  wavenumbers  Sq  upon  wavenumbers  in  the  range  (k-q,  k),  without  loss 
from  the  wavenumbers  <q. 

The  most  striking  feature  of  (6.3)  is  that  is  about  eight  times  as 
large  as  v^.  Some  help  in  interpretation  comes  from  noting  that  -  v,_ 
arises  from  the  first  term  on  the  right-hand  side  of  (2.7),  or  (2.8). 

According  to  the  discussion  in  Sec.  4,  this  term  would  contribute  all  of 
v£  if  the  Eulerian  field  were  frozen.  The  remaining  contribution  repre¬ 
sents  the  effects  of  change  of  the  Eulerian  field  with  time.  A  large  ratio  of 
t0  v. 

v' ^means  that  the  migration  of  particles,  under  the  action  of  wavenumbers 
k  * 

>k,  is  accompanied  by  accelerations  which  leave  the  large-scale  Eulerian 
field  almost  unchanged;  the  particles  are  subjected  to  pressure  forces 
such  that  they  take  on  nearly  the  velocity  of  their  new  surroundings,  apart 
from  fluctuations  at  scales  Sk  }  In  other  words,  the  transport  of  particles 
by  the  velocity  excitation  in  wavenumbers  >k  is  much  more  efficient  than 
the  transport  of  momentum  across  gradients  of  the  large-scale  velocity  field. 

The  migration  of  fluid  particles  is  directly  related  to  the  turbulent 
transport  of  heat  in  a  fluid  with  passive  temperature  gradients.  A  large 
ratio  v£/\>k  therefore  implies  that  the  ratio  of  eddy  conductivity  to  eddy 
viscosity  is  large.  The  LHDI  result  for  the  asymptotic  thermometric  eddy 
conductivity  due  to  wavenumbers  >k  is,  in  fact,  just  v£  -  v^.  Thus  (6.3) 
implies  that  the  eddy-transport  Prandtl  number  is 


°k  =  V(vk  -  V  =  °-14 


(inertial-range  k). 


(6.4) 


Both  G(t)  and  R(t)  have  zero  slopes  at  t  =  0.  The  result  for  G  is 
obvious  from  (5.4).  The  slope  of  R  vanishes  because  the  right-hand  side 


[ 
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of  (5.5)  is  proportional  to  the  energy  transfer  if  t  =  0  [see  (4.2)],  and 
the  latter  vanishes  in  the  inertial  range.  Figure  2  indicates  that  G  and 
R  are  indistinguishable  near  the  origin  but  are  markedly  different  at 
large  argument  values.  For  t  <,  1,  the  numerical  investigation  supports 
the  analytical  result  that  the  behavior  of  G(x)  and  R(t)  is  dominated  by 
interactions  local  in  wavenumber.  The  behavior  at  large  t  turns  out  to 
be  dominated  by  interactions  with  low  wavenumbers.  Analysis  of  this  phe¬ 
nomenon  requires  investigating  the  asymptotic  forms  of  all  the  terms  on 
the  right-hand  sides  of  (5.4)  and  (5.5).  Only  the  results  and  their  physical 
interpretation  will  be  described  here. 

Consider  (5.5)  first.  The  low-wavenumber  contributions  come  from 
small  v  and  small  w.  The  small-w  contributions  turn  out  to  be  negligible 

for  all  t.  The  dominant  small-v  contributions  for  large  t  come  from 

-3/2  .  . 

v  i  t  and  appear  in  the  second  and  sixth  terms  on  the  right-hand  side 

of  (5.5),  both  of  which  have  the  coefficient  Dyw.  These  contributions 

act  as  a  source  term  which  supports  R(x)  against  the  effective  eddy  damping 

from  v,  w  =  0(1).  The  result  is 

R(x)  -v  x-2  (t  >>  1),  (6.5) 

apart  from  a  factor  which  varies  more  slowly  than  any  power  of  x.  In  (5.4), 

-3/2 

the  important  contributions  again  come  from  v  >,  x  and  appear  in  the 
second  and  last  terms  on  the  right-hand  side.  They  result  in  a  tail 

G(x)  -v  x"15/2  (t  »  i),  (6.6) 

apart  from  a  factor  which  varies  more  slowly  than  any  power  of  x.  The 

numerical  results.°are  consistent  with  (6.5)  and  (6.6).  The  dominant  con- 
A 
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tributions  to  C,  as  given  by  (5.11),  and  to  the  integrals  in  (6.2),  are  all 
from  values  of  the  time-arguments  such  that  G  and  R  are  dominated  by  inter¬ 
actions  local  in  wavenumber. 

The  physical  interpretation  of  the  tails  on  R  and  G  is  interesting. 
Transfer  of  energy  from  low  to  high  wavenumbers  involves  straining  pro¬ 
cesses  which  produce  relative  displacements  of  fluid  elements.  Suppose 
that  there  is  an  initial  excitation  of  the  velocity  field  at  wavenumber  p 
and  that  this  excitation  is  spatially  distorted  by  the  action  of  much 
higher  wavenumbers  at  a  time  t^.  Apart  from  the  effects  of  acceleration, 
spectrum  changes  can  be  produced  just  by  the  differential  displacement. 

If  the  elements  along  some  line  in  the  fluid  suffer  displacement  modulated 
at  wavenumbers  'v-k  (k  >>  p),  then  the  resulting  Eulerian  velocity  field 
will  have  a  component  at  wavenumber  k  because  elements  which  are  equally 
spaced  after  the  displacement  were  not  equally  spaced  before  displacement. 

Of  course,  there  will  be  additional  spectrum  changes  because  the  fluid 
elements  suffer  changes  of  velocity,  but  they  turn  out  not  to  be  relevant. 

Now  consider  u(x,t|r)  for  t  >  t^  and  r  <  t^.  The  field  u(,x,t|t)  is 
the  Eulerian  field  after  distortion,  and  has  a  component  k.  The  field 
u(x,t|r)  is  the  Eulerian  field  before  distortion,  but  relabeled  according 
to  the  positions  of  the  fluid  elements  after  the  distortion.  Because  of 
the  modulated  displacement,  correlation  in  the  velocity  of  a  fluid  element 
before  and  after  distortion  shows  up  as  a  correlation  between  the  fields 
u(x,t|t)  and  u(x,t|r)  that  contributes  to  U(k;t|r).  But  for  r  <  t^,  the 
excitation  of  interest  is,  in  the  Eulerian  picture,  at  the  low  wavenumber  p. 


The  connection  with  wavenumber  k  comes  only  through  the  labeling  transforma- 


-27- 


tion  that  takes  place  between  t^  and  t.  Consequently,  there  is  a  k  component 
of  u(x,t|r)  that  varies  with  r  at  the  slow  rates  of  evolution  of  the  low 
wavenumber  p,  and  the  contribution  of  the  distortion  process  to  U(k;t|r) 
decays  slowly  with  r. 

The  terms  in  (5.5)  which  give  rise  to  the  tail  in  R  can  be  seen  to 
describe  just  this  process  if  the  perturbation  analysis  is  traced  through 
in  detail.  The  contribution  to  U(k;t|r)  depends  on  the  shear  associated 
with  the  low  wavenumbers  p;  unless  velocity  gradients  are  present,  the 
modulated  displacement  of  fluid  elements  cannot  generate  modulation  of 
the  velocity  field. 

A  similar  but  more  involved  argument  gives  the  physical  basis  for  the 
tail  on  G.  The  essential  point  is  that,  because  of  differential  displace¬ 
ment  of  fluid  elements  (labeling  transformation),  the  externally  imposed 
perturbation  in  the  field  u(x,t|r)  which  produces  a  pure  wavenumber  k 
perturbation  in  u(x,t|t)  must  contain  all  wavenumbers,  including  low  wave- 
numbers  p  which  evolve  slowly. 

According  to  the  LHDI  closure  formulas,  the  functions  R(k;t|r)  and 
G(k;t|r)  determine  the  effective  memory  times  for  the  buildup  of  triple 
correlations  and  energy  transfer.  The  action  of  wavenumbers  <<k  in  making 
the  functions  more  persistent  therefore  tends  to  lengthen  these  memory 
times.  An  interpretation  is  the  following:  Production  of  high  wavenumbers 
involves  stretching  of  large-scale  spatial  structures  into  structures  having 
small  transverse  scales.  During  this  process,  the  correlations  present 
in  the  large  scales  are  not  destroyed  exponentially  wii:h  repeated  halvings 
f  scale  size,  but  tend  to  persist  throughout  the  straining  process.  In 
wavenumber  language,  the  triple  correlations  associated  with  energy  transfer 
tend  themselves  to  be  transferred  up  the  wavenumber  spectrum. 
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7.  LAGRANGIAN  AND  EULERIAN  STRUCTURE  FUNCTIONS 

A  generalized  spacetime  structure-function  tensor  may  be  defined  by 

B.j(x,t|r;x',t' |r')  =  ^[u.(x,t|r)  -  u^x*  ,t’  |r' )] 

x  [Uj (x,t |r)  -  u_.(x' ,t' |r')]).  (7.1) 

The  structure  tensor  measures  the  velocity  difference  between  two  fluid 
elements.  The  invariance  property 

3U„  (x ,t  | r ;x,t  | r’  )/3t  =  0  (7.2) 

for  homogeneous  turbulence  was  obtained  in  Ref.  1.  It  expresses  that  the 
transformation  from  Eulerian  to  Lagrangian  coordinates  preserves  volume  in 
an  incompressible  fluid.11  For  isotropic  turbulence,  (7.1)  and  (7.2)  give 

Bij(x,t|r;x’,t'  |r')  =  {.^(r)  +  U(r')]  -  2U.  j  (x.t  |r;x'  ,t’  |r' ) ,  (7.3) 

where  U(r)  is  the  mean-square  velocity  in  any  direction  at  time  r. 

A  purely  Lagrangian  structure  tensor,  which  compares  fluid  elements 
labeled  simultaneously,  is  given  by  t  =  t’,  and  a  purely  Eulerian  tensor 
is  given  by  r  =  t,  r'  =  t'.  The  specialized  Lagrangian  tensor 

B^-j (*,-£,'  >t  I r )  =  B„(x,t|t;x'  ,t|r)  (7.4) 

can  be  evaluated  for  t  >  r  from  the  abridged  LHDI  equations.  This  tensor  is 

purely  solenoidal,  while  the  full  Lagrangian  structure  tensor  is  not. 

When  the  turbulence  is  isotropic  and  stationary,  standard  Fourier 
12 

transformation  formulas  yield 

B(y;t|r)  =  4/  [1  -  S'1^(k^)  R(k;t  |r)]E(k)dk, 

0  y 
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B|.(y;t|r)  =  4/  {-i  -  [Sin(ky)_  _  cos(ky)  ]R(k.t  |r)}E(k)dkt 
0  3  (ky)3  (ky)Z 


where  B  and  B|(  are,  respectively,  the  trace  of  B„(y;t|r)  and  the  diagonal 
component  parallel  to  y,  in  a  coordinate  system  aligned  with  y.  The  diag¬ 
onal  components  perpendicular  to  y  are  each  given  by  Bx  =  (B  -  B||)/2,  and 
the  off-diagonal  components  vanish  in  this  coordinate  system. 

If  the  inertial-range  forms  (5.1)  and  (5.2)  are  substituted  into  (7.5), 
the  result,  denoted  here  by  superscript  (1),  is 

B<l)(„t|r>  =  t2/3y2/3F<l,<t1/3T/y2/3), 

B,!1)(yit|r)  =  .2/3y  2/3r!1)(.1/3T/y2/3),  ,  .  t  -  r.  (7.6) 


where 


F(1)(s)  =  4C /  [1  -  SiEJi  R(x2/3s)]x‘5/3dx, 
0  X 


F,  ;(s)  =  4C/  [i  -  ( 


1  .sin  x  cos  x  .D,  2/3  -5/3 


)R(x  '  s) ]x  dx. 


and,  in  particular, 

F(1)(0)  =  |  r(-~)  C  =  8.52, 


F^1)(°)  =  f(1>(0).  (7.7a) 


The  computed  value  of  C  is  used  in  (7.7a). 

The  integrations  over  x  in  (7.7)  are  transformed  wavenumber  integrations, 
and  they  converge  at  both  limits  as  a  consequence  of  the  properties  of  R 
established  in  Secs.  5  and  6.  And  if  y  is  the  reciprocal  of  a  wavenumber 
within  an  extensive  inertial  range,  only  inertial-range  wavenumbers  contri¬ 
bute  appreciably  to  the  integrals.  Nevertheless,  (7.6)  does  not  give  the 
inertial-range  structure  functions  correctly.  This  can  be  seen  from  (4.2). 

If  (4.2)  is  integrated  over  the  energy-containing  range,  it  yields 
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/  [3R(k;t|r)/3tl  E(k)dk  =  -e/2.  (7.8). 

4  1  t=r 

energy 

range 

Consequently,  if  t  =  t-r  is  in  the  inertial  range  of  time  differences, 
small  compared  to  the  characteristic  evolution  times  of  the  energy-contain¬ 
ing  wavenumbers. 


/  R(k;t|r)E(k)dk  =  -  (t-r)e/2.  (7.9) 

energy 
range 

The  factor  ky  in  (7.5)  can  be  replaced  by  unity  with  negligible  error  if 
k  is  in  the  energy-containing  range  and  y  ^  in  the  inertial  range.  Hence 
there  is  an  additional  contribution 

B(2)(y;t|r)  =  2(t  -  r)e,  2 J ( y  ;t|r)  =  2(t  -  r)e/3,  (7.10) 


so  that 


B(y jt | r)  .  e2/V/3«el/3t/y2/3),  B„(y;,|rt  .  ^V'\u1/3,/y2/h , 


with 


F(s)  =  F(1)(s)  +  2s, 


t  =  t  -  r  >  0, 


F,  (s)  =  F(1)(s)  +  2s/3. 


(7.11) 

(7.12) 


If  s  ^  1,  the  two  contributions  to  F(s)  have  comparable  magnitudes. 

There  is  no  corresponding  asymptotic  contribution  from  dissipation- 
range  wavenumbers.  The  total  transfer  into  the  dissipation  range  equals 
the  transfer  out  of  the  energy  range,  but  the  dissipation-range  contribution 
to  (7.5)  is  suppressed  by  the  trigonometric  coefficients  if  y  is  in  the 
inertial  range  and  by  the  small  correlation  time  of  R(k;t|r)  if  t-r  is  in 
the  inertial  range. 

A  closer  look  at  how  the  contributions  (7.10)  arise  shows  that  actually 
they  should  not  be  considered  direct  effects  of  the  energy-containing  eddies. 
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Equation  (7.2)  can  be  Fourier  analyzed  just  like  the  energy-balance  equation. 
It  then  describes  a  similar  cascade  process:  There  is  a  negative  contribu¬ 
tion  -26„e/3  to  the  left-hand  side  from  the  energy-containing  range,  if 
t-r  is  small  compared  to  the  characteristic  times  of  the  energy-containing 
eddies.  This  is  balanced  by  a  positive  contribution  at  higher  wavenumbers. 
These  properties  can  be  established  by  analysis  similar  to  that  which  g.ves 
(4.2).  Note  also  the  discussion  of  small  scales  at  the  end  of  Sec.  1. 

Thus,  if  (7.1)  is  Fourier  analyzed  before  (7.2)  is  used,  there  are  two 
equal  and  opposite  energy-range  contributions  which  cancel,  and  (7.10) 
appears  as  a  contribution  at  higher  wavenumbers,  from  the  cascade  process. 

The  formal  analysis  evidently  indicates  an  ambiguous  origin  for  (7.10), 
determined  by  what  implicit  uses  of  conservation  and  invariance  properties 
have  been  made. 

The  physical  interpretation  of  (7.10)  seems  fairly  straightforward: 

On  the  average,  two  particles  which  are  separated  by  y  at  a  given  time  t 

are  undergoing  a  relative  acceleration  negatively  correlated  with  the 

velocity  difference,  the  correlation  being  measured  by  the  slope  of 
/  \ 

Bv  ;(y;t|r)  as  a  function  of  r  at  r  =  t.  This  was  discussed  at  the  end  of 
Sec.  1.  ^  (y ;t | r)  has  zero  slope  at  r  =  t  and  is  a  contribution  from  the 

part  of  the  relative  acceleration  which  is  uncorrelated  with  the  relative 
velocity  at  time  t.  The  correlation  described  by  B'  '(y;t|r)  is  due  to 
straining  processes  established  over  finite  times.  If  the  velocity  is 
multivariate  Gaussian  at  tQ,  so  that  T(k,tQ)  vanishes  for  all  k,  then  there 
is  no  initial  correlation  of  relative  acceleration,  due  to  pressure  forces, 
with  relative  velocity. 

The  abridged  LHDI  equations  of  the  present  paper  do  not  predict 
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B(y;t|r)  for  r  >  t.  But  clearly  the  negative  slope  as  a  function  of  r 
cannot  persist.  Eventually,  B(y;t|r)  must  increase  with  r.  Note  that 
if  two  particles, known  to  be  separated  by  y  at  t,  have  reached  a  new 
separation  y'  at  a  later  time  t',  their  mean-square  velocity  difference  is 
not  given  by  B(y * ;t  * 1 1  * ) .  These  two  particles  are  members  of  a  special 
subset  of  all  pairs  separated  by  y'  at  t',  and  have  special  statistics. 

The  present  results  at  r  =  t  thus  do  not  contradict  the  proposal  by  C.  C. 

Lin  that  the  simultaneous  correlation  of  relative  velocity  and  relative 
acceleration  be  positive  and  constant  at  times  sufficiently  long  after 
release  of  a  particle  pair  at  an  inertial-range  spatial  separation.  The 
structure  function  which  describes  the  change  with  time  of  simultaneous 
relative  velocity  is  B .  ^x.t  1 1 '  ;x-y ,t  1 1 ' ) . 

Equation  (7.11)  is  inconvenient  for  finding 

B(t-r)  i  B(0;t|r). 

which  gives  the  mean-square  change  of  velocity  of  a  single  fluid  element 

in  the  interval  t-r.  For  inertial-range  values  of  t  (not  for  t  0), 

00 

B(t)  =  ke|t|,  k  =  6Cj  [1  -  R(s)]s  ^ds  +  2.  (7.13) 

0 

The  second  term  in  tc  is  obtained  from  (7.10),  and  the  first  term  is  obtained 
by  using  (5.1)  and  (5.2)  v  (7.5)  with  y  =  0.  Equation  (7.12)  yields 

F(s)  as  ks  F  (s)  *  ks/3,  s  »  1.  (7.12a) 

B(t)  must  be  even  in  t  for  stationary  turbulence,  and  this  fact  has  been 

incorporated  into  (7.13).  The  linear  dependence  of  B(t)  on  t  in  the  inertial 

14  ... 

range  was  obtained  by  Inoue  by  Kolmogorovian  dimensional  analysis. 

The  computed  results  for  F(s)  and  ^(s)  are  plotted  in  Fig.  4.  The 


! 
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value  computed  for  the  integral  in  (7.13)  is  1.46,  yielding  k  =  17.5. 

The  above  formulas  for  the  inertial-range  Lagrangian  structure  func¬ 
tions  assume  stationary  isotropic  turbulence,  which  must  be  maintained  by 
statistically  stationary  driving  forces  of  some  kind.  The  results  are  un¬ 
changed  if  the  turbulence  is  freely  decaying.  This  can  be  verified  by  re¬ 
tracing  the  analysis  carefully  for  the  freely  decaying  case.  It  also  can 
be  seen  by  a  simple  argument:  The  decay  of  the  total  kinetic  energy  is 

O  , 

given  by  3u  /3t  =  -2e,  where  u  is  the  total  rms  velocity.  The  characteristic 

2  2 

decay  time  of  u  is  therefore  u  /e.  In  a  time  t  <<  u  /e,  the  mean-square 
change  of  velocity  due  to  decay  (i.e..  due  to  the  lack  of  energy-range 
sustaining  forces)  is  ({Au)  (et/u)  et[t/(u  /e)],  which  is  asymptotically 
negligible  compared  to  the  mean-square  velocity  differences  of  order  et 
found  above.  Note  that  mean-square  change  of  velocity  does  not  mean  change 
of  mean-square  velocity. 

A  purely  Eulerian  structure  tensor  is  obtained  by  setting  r  =  t, 

r'  =  t'  in  (7.1).  For  t  -  t'  small  enough  that  the  time  dependence  of 

E 

E(k,t)  is  negligible,  the  formulas  for  the  Eulerian  scalars  B  (y;t,t') 
and  BE(y;t,t' ),  analogous  to  B  and  B|( , are  identical  with  (7.5),  except  that 
R(k;t|r)  is  replaced  by  the  Eulerian  modal  correlation  function 

RE(k;t,t’ )  =  27rk2U(k}t|t;t'|t')/[E(k,t)E(k,t')]1/2.  (7.14) 

The  Lagrangian  and  Eulerian  modal  time  correlations  measure  quite 
different  properties  when  k  is  in  tue  inertial  range.  The  Lagrangian  cor¬ 
relation  is  determined  by  intrinsic  distortions  of  the  small-scale  structures, 
and  the  correlation  time  is  of  order  (v^k)”1,  where  vk  is  the  rms  velocity 
in  wavenumbers  >  k.  The  Eulerian  correlation  is  dominated  by  the  convection 
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of  the  small  scales  by  the  energy-containing  scales. 

The  convection  of  a  frozen  small-scale  structure  by  a  uniform  velocity 
u  induces  the  time  dependence  exp(-ik*u)  in  the  Fourier  amplitude  of  the 
small-scale  structure.  Now  assume,  in  accord  with  Kolmogorov's  hypotheses, 
that  the  distribution  of  the  small-scale  velocity  fluctuations  is  statisti¬ 
cally  independent  of  the  simultaneously  measured  distribution  of  the  con- 
vecting  velocity.  Then  it  follows  immediately  that 

RE(k;t,t--r)  =  M(kr,t),  (7.15) 

where 

M(a,t)  =  (e-12-’-)  =  /e_ia#J,P(u,t)d3u  (7.16) 

is  the  characteristic  function  of  the  one-point  probability  distribution 
P(u,t)  of  the  convecting  velocity  at  time  t.  The  correlation  time  indicated 
by  (7.15)  is  (vQk)  \  where  vQ  is  the  rms  value  of  any  vector  component  of 
the  convecting  velocity,  now  identified  with  the  energy-range  velocity. 

Since  this  correlation  time  is  <<  (v^k)-1,  it  may  be  concluded  that  dis¬ 
tortion  effects  make  a  negligible  change  in  the  Eulerian  time  correlation 
and  that  (7,15)  is  the  complete  asymptotic  result  for  inertial-range  k. 
Finally,  P(u,t)  in  (7.16)  may  be  taken  as  the  one-point  distribution  of 
the  entire  velocity,  since  the  energy-range  contributions  dominate  the 
latter  distribution. 

Equation  (7.15)  follows  also  for  dissipation-range  wavenumbers  if 

the  relevant  distortion  time  for  scales  of  order  k  E  is  assumed  to  be  not 

1  '2 

less  than  the  Kolmogorov  time  (v/e)  '  regardless  of  how  high  k  may  be.  In 
other  words,  it  is  necessary  that  the  very  small-scale  structures  follow 
the  structures  of  scale  kg  ^  and  not  exhibit  independent  distortions. 
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The  result  for  the  Eulerian  modal  response  function  G  (k;t,t') 
corresponding  to  (7.15)  is 

GE(k;t,t-T)  =  exp( -vk2T)M(kr  ,t ) ,  (7.17) 

for  k  in  the  inertial  or  dissipation  range.  The  viscous  factor  appears  in 
(7.17)  because  an  arbitrary  perturbation  at  a  wavenumber  high  in  the  dissi¬ 
pation  range  will  be  uncorrelated  with  the  existing  structures,  will  not 

be  supported  by  them,  and  will  decay,  during  the  convection  by  energy-range 

2 

scales,  as  if  in  isolation.  Note  that  vk  <<  v.k  ,  so  that  when  the  mer- 

s  u  s 

tial  range  is  extensive  the  viscous  damping  in  (7.17)  is  negligible  except 
for  k  in  the  far  dissipation  range. 

Equations  (7.15)  and  (7.17)  are  presumably  exact  results  if  the 

Gil 

underlying  assumption  of  asymptotic  statistical  indepdence  of  simultaneously 
measured  small-  and  large-scale  structures  is  valid.  No  closure  approxi¬ 
mation  is  used,  and  isotropy  of  the  energy-containing  range  is  not  invoked. 
The  abridged  LHDI  equations  do  not  predict  Eulerian  time  correlations. 

However,  the  unabridged  LHDI  approximation  leads  easily^  to  explicit  results 

E  E  -I 

for  R  (k;t,t')  and  G  (k;t,t')  when  k  is  high  enough  that  (vQk)  is  short 

compared  to  Lagrangian  correlation  times  for  mode  k,.  In  that  case,  the 

Galilean  invariance  properties  of  the  LHDI  equations^  imply  directly  that 

RE(k;t,t-t)  =  exp[-  jV^kV], 

GE(k;t,t-x)  =  exp[-vk2t  -  j  v2k2t2],  (7.18) 

in  the  isotropic  case.  This  agrees  with  (7.15)  and  (7.17)  if  the  univariate 
velocity  distribution  is  Gaussian.  Thus  the  accuracy  of  the  LHDI  results 
depends  on  how  close  the  actual  univariate  distribution  is  to  Gaussian. 


passra*?***1- 


i; 

i 
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In  the  case  of  inhomogeneous  turbulence  with  nonzero  mean  velocity,  the 
LHDI  result  again  is  exact  if  and  only  if  the  univariate  velocity  distri¬ 
bution  is  normal. 

E  E 

Insertion  of  (7.18)  into  the  formulas  for  B  and  yields 
BE(y;t,t* )  =  e2/3y2/3FE(v  T/y),  BE(y;t,t')  =  e2/3y2/3FE  (v.-r/y) , 


T  =  t  -t'. 


where 


„E ,  ,  r ,  sin  x  ,12  2,,  -5/3, 

F  (s)  =  4Cj  [1 - exp(-  x  s  )]x  dx, 

«  x  z. 


(7.19) 


F,f(s)  =  4C/  [  \  -  (^ 


sin  x  cos  x  ,  ,  12  2,-,  -5/3, 

g — - 2 —  >  exp(-  -  x  s  )]x  dx. 


(7.20) 


Also, 


BE(t)  =  xe(ev  t)2/3,  tce  =  4C /  [1  -  exp(-  i-  s2)]s_5/3ds. 


(7.21) 


E  E 

where  B  (t-t')  =  B  (0;t,t').  Equation  (7.21)  gives  the  mean-square  velocity 

change  at  a  fixed  point  in  laboratory  coordinates.  Equations  (7.19)-(7.21) 

are  valid  for  y  and  vqt  in  the  inertial  range  of  spatial  scales.  There  are 

£ 

no  contributions  analogous  to  (7,10);  R  (k;t,t')  has  zero  slope  at  t  =  t' 

for  all  wavenumbers,  in  stationary  turbulence. 

E  E 

F  (s)  and  F)(  (s)  are  plotted  in  Fig.  5.  The  integral  in  (7.21)  has  the 
value  1.61,  yielding  <e  =  11.4.  Equation  (7.20)  yields 

FE(s)  w  <es2/3,  FE(s)  w  Kes2/3/3,  s  »  1. 


(7.20a) 
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8.  LAGRANGIAN  ACCELERATION  AND  PRESSURE  COVARIANCES 
The  function 

a(x,t|r)  =  3u(x,t|r)/3r  (8.1) 

is  the  acceleration  at  time  r  of  the  fluid  element  which  arrives  at  :x  at 
time  t.  The  two  simplest  Lagrangian  covariances  between  acceleration  and 
velocity  are 

K_(y;t|r)  =  (u.  (x,t  |t)a..  (x-y  ,t \rp  (8.2) 

and 

K?j(y;t|r)  =  ^(x.t  |r)a..  (x-y  ,t  |t)^.  (8.3) 

Giver,  a  pair  of  particles  known  to  be  separated  by  y  at  time  t,  K„  measures 
the  correlation  between  the  velocity  of  one  at  time  t  and  the  acceleration 
of  the  other  at  a  time  r.  measures  the  correlation  between  the  vel¬ 

ocity  of  the  first  at  time  r  and  the  acceleration  of  the  second  at  time  t. 

The  relation  between  K._.  and  kT..  is  simplest  in  stationary,  homogeneous 
turbulence.  To  permit  such  turbulence,  let  a  stationary,  homogeneous, 
solenoidal  forcing  field  f(x,t)  be  introduced  so  that  the  Navier-Stokes 
equation  is 

[a /at  +  u(x,t jt)-V]u(x,t|t)  =  a(x,t|t),  (8.4) 

a(x,t | t)  =  vV2u(x,t|t)  -  Vp(x,t)  +  f(x,t),  (8*5) 

where  p(x,t)  is  the  kinematic  pressure  field.  The  nature  of  f(x,t)  need 
not  be  specified  further.  The  forces  may,  but  need  not,  be  negative  damp¬ 
ings  as  in  Sec.  4.  In  the  stationary  state. 


K. . (x-x* ;t|r)  =  3U. . (x,t | t ;x' ,t |r)/3r  =  -  3U. .(x,t |t;x' ,t |r)/3t.  (8.6) 
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An  expression  for  the  last  member  of  (8.6)  is  readily  obtained  by  multi¬ 
plying  (8.4)  with  u(x',t|r),  multiplying  (1.2)  with  u(x,t|t),  adding,  and 
averaging.  This  gives 

Ki j  <  £-x' ;t | r)  =  (uj(x',t|r)un(x,t|t)3u.(x,t|t)/3xn> 

+  ^ui(x,t|t)un(x' ,t|t)3Uj(x' ,t|r)/3x^)  -  ^u.(x’ ,t|r)a.(x,t|t)X  (8.7) 

If  the  second  term  on  the  right-hand  side  is  rewritten  using 

3u  (x'.tltJ/Sx1  =  0  and  homogeneity,  (8.7)  becomes 
n  ~  n 

K. ,(x-x' it | r)  a  -  Kt.(x'-x;t|r) 
ij  ~  ]i  ~  — 

+  ^u.(x) ,t|r)[un(x,t|t)  -  un(x' ,t|t)]3ui(x,t|t)/3xn\  (8.8) 

Equation  (8.8)  shows  that  K„(x-x'  ;t|r)  =  - 
but  not  if  x  i  x1 .  The  result  for  x  =  x*  can  also  b'  obtained  from  dif¬ 
ferentiation  of 


U. .(x,t | t ;x,t | r)  =  U. ,(x,r|t;x,r|r),  (8.9) 

J-j  **  ~  ~  - 

a  relation  which  follows  from  (7.2).  Moreover,  the  evenness  of 

U. . (x,t | t;x,t | r)  in  t-r  shows  that  K(t-r)  e  K..(0;t|r)  is  an  odd  function 

of  t-r. 


is  solenoidal  in  i  because  the  Eulerian  velocity  is  solenoidal, 
and  is  solenoidal  in  j  also  if  the  turbulence  is  homogeneous  and  reflection- 
invariant.  In  the  stationary  state, 


-  3U„(x,t  |t;x'  ,t  |r)/3t  =  j  (x,t  1 1  ;x'  ,t  |r)/3t , 

so  that  for  stationary,  isotropic  turbulence  (7.11)  yields 

K(y;t|r)  =  j  eF( e1/3i/y2/3) ,  K„(y;t|r)  =  |  eF(|  ( e1/3i/y2/3) , 


(8.10) 


t  =  t  -  r  >  0 


(stationary  turbulence) 


(8.11) 


-49- 


where  K  and  K(|  are  the  trace  and  the  longitudinal,  diagonal  element  of 

K...  The  functions 
il 

F(s)  =  dF(s)/ds,  F|((s)  =  dF(|  (s)/ds 
are  plotted  in  Fig.  6.  Equation  (7.13)  yields 

K(t)  =  y  <z  sgn  x  (stationary  turbulence).  (8.12) 

It  must  be  remembered  that  (8.12)  is  established  only  for  |t|  in  the  inertial 
range  of  times  and  is  invalid  for  |t|  -*■0. 

The  change  of  sign  in  (8.12)  has  a  simple  physical  interpretation. 

For  t  >  0,  K(t)  measures  the  correlation  of  a  particle's  velocity  with  its 
previous  acceleration.  This  is  positive  because  the  velocity  increment  pro¬ 
duced  by  the  acceleration  shows  some  persistence.  On  the  other  hand,  a 
negative  correlation  of  velocity  with  later  acceleration  is  necessary  if 
the  velocity  is  to  change  direction  without  increasing  its  magnitude,  in 
mean  square.  Contrast  this  behavior  of  one-particle  velocity  and  accelera¬ 
tion  with  the  opposite  correlation  of  relative  velocity  and  relative  acceler¬ 
ation  of  two  particles,  as  indicated  by  the  slope  of  B(y,t|r)  at  r  =  t 
[Equations  (7.10)  and  (7.11),  et  seq.]. 

A  constant  negative  value  of  K(x)  for  t  negative  and  in  the  inertial 
range  has  been  proposed  on  the  basis  of  a  kind  of  Brownian-motion  model  of 
the  behavior  in  velocity  space . ^ ^  The  present  result  does  not  seem 
properly  explicable  in  terms  of  such  a  model.  The  major  contribution  to 
K(t)  at  an  inertial-range  t  comes  from  eddies  whose  circulation  times  are 
of  order  t  .  A  Brownian-motion  model  requires  that  instead  the  major 
effect  come  from  eddies  whose  circulation  times  are  <<  t,  so  that  many  un¬ 
correlated  eddies  contribute.  The  similarity  of  the  results  of  the  two 
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approaches  seems  fortuitous,  but  this  may  be  a  premature  conclusion. 
Conservation  of  energy  requires 


(u^x.tlt^ix.tlt))  =  j  (8.13) 

in  stationary,  isotropic  turbulence.  If  the  driving  forces  f^  are  confined 
to  the  energy-containing  range,  then  (8.13)  is  nearly  satisfied  with 
u^(x,t|t)  replaced  by  u^(x_',t|r),  provided  that  |x-x'  |  and  t-r  are  small 
compared  to  energgy-range  space  and  time  scales.  In  this  case,  (8.5)  yields 


where 


Kt .  (x-x»  Jt:  |  r)  =  W2  U  (:.,t|r;x',tlt)  +  D.  .(x-x’  ;t  |r)  +  i  6..e,  (8.14) 

j_j  X  1 J  -v-  **  1J  **  °  1J 

Dij(x-x'Jt|r)  =  -<u.(x,t|r)3p(x',t)/3x^> 


(8.15) 


is  the  covariance  of  pressure  gradient  with  Lagrangian  velocity.  D„  is  a 
pure  antisolenoidal  tensor  in  homogeneous,  reflection-invariant  turbulence. 
Its  curl  with  respect  to  each  index  vanishes. 

The  right-hand  side  of  (8.15)  can  be  converted  into  a  third-order 
moment  of  the  velocity  field  by  eliminating  the  pressure  in  standard  fashion. 
The  triple  moment  can  be  evaluated  for  t  >  r  by  the  abridged  LHDI  approxi¬ 
mation  of  Ref.  1.  The  result  for  isotropic  turbulence  is 


D(y;t|r)  =  /  Q(k;t|r)sA.ffiI  dk, 

0  *y 

DCy;t|r)  =  J°°Q(k;t (8.16) 
0  (kyr  (ky)2 


where 

Q(k;t|r) 


-4A2JJAdpdq  P2qhkpq 


t 

U(p ;t | r ) /  U(q;t |s)ds, 
r 


h, 

kpq 


=  (l-z2)(z+xy). 


(8.17) 
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In  (8.16),  D  and  Dx  are  the  trace  and  diagonal  element  perpendicular  to^y. 
The  diagonal  element  parallel  to  y  is  given  by  D  =  D  -  2D  ,  and  the  off- 
diagonal  elements  vanish  in  a  coordinate  system  aligned  with  y.  In  (8.17), 
x,  y,  and  z  are  the  cosines  of  the  interior  angles  opposite  k,  p,  and  q  in 
a  triangle  with  these  wavenumbers  as  sides. 

For  |x-x' |  and/or  t-r  in  the  inertial  range  of  separations,  assume 
that  the  v  term  in  (8.14)  is  negligible,  subject  to  later  check.  Then 
is  purely  curlfree,  in  contrast  to  K„  which  is  purely  solenoidal.  If  the 
inertial-range  forms  (5.1)  and  (5.2)  are  substituted  into  (8.17),  the 
result  can  be  written 

Q(k;t | t-T )  =  -C2ek-1Q(e1/3k2/3T),  (8.18) 

where 

Q(s)  =  /  duj  dv  h  u'5/V8/3R(u2/3s)J  R(v2/3w)dw.  (8.19) 

0  |u— 1|  iUV  0 

The  function  Q(s)  is  plotted  in  Fig.  7.  For  s  >>  1,  Q(s)  *  s  ,  to  within 
a  factor  that  varies  more  slowly  than  any  power  of  s.  Equation  (8.18)  gives 

D(y;t|r)  =  -  cZ(c1/3j/y2/3) ,  Dx(y;t|r)  =  -  tZx( e1/3i/y2/3) , 

t  =  t  -  r,  (8.20) 


Also,  „ 

D( 0 ;t  | r)  =  -  |  C2eJ  Q(s)^-  . 

2  0  s 


(8.22) 
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Using  these  values  for  the  pressure  contribution  to  in  the  inertial 

range,  it  can  be  verified  that  the  viscous  contribution  to  (8.14)  is  negli- 

1/2 

gible  in  comparison  if  t-r  is  large  compared  to  the  Kolmogorov  time  (v/e)  . 

The  check  consists  of  assuming  a  long  inertial  range  and  evaluating  the 

total  viscous  contribution  from  the  whole  range  of  wavenumbers  up  to  k  = 

3  1/4 

(e/v  )  ,  using  the  forms  (5.1)  and  (5.2)  and  using  the  properties  of 

R(t)  for  long  t  discussed  in  Sec.  6. 

The  final  inertial-range  results  are 

K+(y;t|r)  =  -eZ(e1/3x/y2/3)  +  e,  K*(y;t|r>  =  -eZx( e1/3x/y2/3)  +  e/3, 

(stationary  turbulence)  x  =  t  -  r,  (8.23) 

00 

K^(t)  =  -  2.  c2ef  Q(s)—  +  e  (stationary  turbulence),  (8.24 
2  0  S 

where  Kty;t|r)  and  K^(y;t|r)  are  the  trace  and  transverse  diagonal  element 

X  XX 

of  K.,(y;t|r),  and  K  (t-r)  =  K  (0;t|r).  The  functions  Z(s)  and  Z  (s)  are 

1 J  *+  * 

plotted  in  Fig.  8.  The  computed  value  of  the  right-hand  side  of  (8.24) 
is  -8.96  e.  Equation  (8.23)  yields 

K+(y;t|t-x)  »  K+(t),  K^(y;t|t-x)  ~  j  K+(t),  e^x/y2^3  >>  1.  (8.23a) 

According  to  (8.8), 

K(x)  =  -  K+(x).  (8.25) 

Thus  if  (7.13),  (8.12),  and  (8,24)  are  consistent  equations,  C  and  R(s) 
must  satisfy 

00  00 

3Cf  [1  -  R(s)  ]s”2ds  +  1  =  |  C2/  Q(s)—  -  1,  (8.26) 

0  2  0  S 

with  Q(s)  given  by  (8.19).  This  consistency  is  assured  by  the  abridged 
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LHDI  procedure.  The  calculation  of  K„  by  (8.6)  amounts  to  evaluating 
the  right-hand  side  of  (8.7),  using  the  abridged  LHDI  recipe  for  all  the 
triple  moments.  The  vanishing  of  tt-e  last  term  on  the  right-hand  side  of 
(8.8)  at  x  =  x'  survives  the  approximation  because  u  (x,t|t)  and  un(x',t|t) 
receive  similar  treatments.  Equation  (8.26)  is  a  useful  overall  check  on 
the  algebra  and  numerical  calculations  (see  Appendix). 

The  Lagrangian  acceleration  covariance 


.  (x-x*  ;t  |  r)  =  ^a^(  x,t  |r)a .  (x' ,t  1 1  ))> 


(8.27) 


i]  ’  '  '  i  1  3 

can  be  evaluated  by  using  the  relation 

A„(x-xj  ;t  |r)  =  3K^(x-x' ;t  |r)/3r. 

This  yields  the  inertial-range  results 

„  |  ,  4/3  -2/34,  1/3  .  2/3,  .  ,  ,  4/3  -2/3-  .  1/3  .  2/: 

A(y;t | t-t )  =  e  y  Z(e  t/y  ),  A  (y;t |t-i)  -Cy  Z  (e  i/y 


(8.28) 


(8.29) 


where  A  and  A  are  the  trace  and  transverse  diagonal  element  of  A...  The 
x  1] 

functions 

Z(s)  =  dZ(s)/ds,  2x(s)  =  dZx(s)/ds 

are  plotted  in  Fig.  9.'  The  quantity  dQ(s)/ds  determined  by  (8.19)  can  be 
evaluated  analytically  at  s  =  0  [see  (8.36),  (8.37)],  when  e 

Z(°)  =Ig[|r(|)]2C2  =  8.39,  Zx(0)=yZ(0).  (8.30) 


The  covariance  of  acceleration  along  a  particle  trajectory  is 
A(t-r)  s  A(0;t|r).  It  satisfies 

A(t)  =  -  dK+(T)/dx.  (8.31) 


Differentiation  of  (8.24)  with  respect  to  t  gives  zero.  This  says  that 
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there  is  no  inertial-range  contribution  to  A(t)  proportional  to  ex-1,  the 
form  that  might  be  suggested  by  Kolmogorov  scaling.  For  given  t,  the 
positive  contributions  to  A(t)  from  some  inertial-range  wavenumbers  are 
balanced  by  negative  contributions  from  other  wavenumbers.  A  relevant  fact 
is  that  the  pressure  forces  conserve  energy  and  cannot  change  the  mean- 
square  particle  velocity.  If  there  were  a  nonzero  net  inertial-range  con¬ 
tribution  of  the  form  ex  \  then  the  inertial  range  (where  only  pressure 

*00 

forces  cause  acceleration)  would  give  a  nonzero  contribution  to  JQA(s)ds, 

the  total  rate  of  change  of  kinetic  energy.  Also,  a  contribution  of  the 

-1  1/2  1/2 
form  et  to  A(x)  would  imply  a  contribution  of  form  -  e  ln(e  x/v  ) 

to  K+(x),  which  would  violate  Kolmogorov's  concept  of  a  v-independent 

inertial  range. 

In  contrast  to  the  results  for  the  structure  functions,  the  inertial- 
range  formulas  for  K„  and  kT^  are  altered  if  the  turbulence  is  freely 
decaying  instead  of  stationary.  The  free-decay  results  are  obtained  simply 
by  subtracting  off  the  driving  term  in  (8.14),  seq.  Thus  (8.11),  (8.12), 
(8.14),  (8.23)-(8,25 )  are  replaced  by 

K(y»t|r)  =  j  eF( z1/Zi/y2/Z)  -e , 

K||(y;t  |  r)  =  ±  eF|(  (  e1/3x/y2/3)  -  e/3,  (8.11') 

K(x)  =  e(-l  t  j  k  sgn  x),  (8.12) 

.  2 

k!.(x-x';t|r)  =  W  U  (x,t|r;x’,t|t)  +  D  (x-x';t|r),  (8.14') 

1 J  ~  ~  X  — 

K+(y;t|t-x)  =  -  eZ(e1,/3x/y2'/3),  K^(y;t|t-x)  =  -  eZx( a1,/3x/y2^3) ,  (8.23') 
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00 

K+(t)  =  -  f  C2e/  Q(s)—  ,  (8.24') 

*  0  s 

K(t)  +  2e  =  -K+(t),  (8.25') 

where  a  corresponding  subtraction  of  the  driving-force  contribution  has  been 
made  from  (8.11)  and  (8.12).  These  changes  express  that  the  decay  of  the 
turbulence  requires  a  deceleration  of  the  particles,  on  the  average.  By 
conservation  of  energy,  the  deceleration  must  contribute  -e  to  K(0),  an 
amount  which  is  comparable  to  the  inertial-range  contributions  and  super¬ 
posed  upon  them.  Equations  (8.18)  -  (8,22)  and  (8.29)-(8.31)  are  unchanged. 

2  2 

The  contribution  of  the  decay  deceleration  to  A.  ^  is  of  order  e  /vq  ,  where 
Vq  is  the  total  rms  velocity  component.  This  is  smaller  than  the  inert ial- 
range  contributions  (8.29)  by  a  factor  of  order(y/iQ)  ,  where  is  the 
characteristic  energy-range  spatial  scale,  and  therefore  does  not  appear 
in  the  asymptotic  formulas. 

It  is  a  matter  of  taste  whether  the  overall  deceleration  in  decaying 
turbulence  should  be  considered  an  energy-range  effect.  Actually,  all 
scales  of  motion  are  involved:  The  loss  of  kinetic  energy  occurs  at 
dissipation-range  scales,  in  the  Eulerian  picture,  so  that  the  entire  cas¬ 
cade  plays  a  role.  It  is  clear,  however,  that  the  scales  of  order  y  play 
no  special  role  in  determining  the  deceleration  contribution  to  K(y;t|r) 
when  y  is  an  inertial-range  separation.  Therefore  it  seems  justified  to 
consider  the  steady-state  formulas  as  the  "true"  inertial-range  expressions. 

The  functions  discussed  above  give  some  limited  information  about 
pressure  statistics.  By  (8.1)  and  (8.15), 
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[3Di^(x-x' ;t|r)/Sr]r=t  =  -  ^ai(x,t|r)3p(x' ,t)/3x^)  . 


(8.32) 


In  homogeneous,  reflection-invariant  turbulence,  only  the  curlfree  term 
-Vp  in  (8.5)  can  contribute  to  (8.32).  Hence, 

^[3p(x,t)/3x^][3p(x|  ,t)/3xj  ])  =  [SD^U-x1  ;t|r)/3r]r_t  .  (8.33) 

The  pressure  spectrum  function  [satisfying  mean-square-pressure  =  j”p(k,t)dk  ] 


is  therefore 


P(k,t)  *  k“^C3Q(kit|r)/3r]r=t 
in  isotropic  turbulence.  By  (8.17),  this  gives 

2 


(8.34) 


P(k,t)  =  //Adpdq  p  qhkpqE(p,t)E(q,t). 


(8.35) 


This  abridged  LHDI  result  for  P(k,t)  is  identical,  after  trigonometric 

manipulations,  with  the  quasinormality  approximation  for  P(k,t)  used  by 
3  12 

Heisenberg  and  Batchelor.  The  unabridged  LHDI  expression  for  P(k,t) 
does  not  appear  to  reduce  to  the  quasinormality  approximation. 


In  the  inertial  range,  (8.35)  yields 


u+1 


P(k,t)  =  C2e4/V7/3J  duj  dv  hl  u-5/V8/3. 

0  |u-l|  iuv 


(8.36) 


17 

Using  trigonometric  identities  and  results  of  Oboukhov,  the  quadratures 
can  be  performed  to  give 


FCk.t)  ■  gg  [  |  r(|>]3C2e4/3k"7/3  .  4.13  «*'V7'3. 


(8.37) 


and 


([p(x,t  -  p(x+y,t)]2,)  =  [Bjj  (jr,t  |t)]2  =  5.39  E4/3y4/3,  (8.38) 

A  Lagrangian  pressure  field  p(x,t|r)  can  be  defined  as  the  pressure 
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at  time  r  in  the  fluid  element  which  arrives  at  x  at  time  t.  The  contribu- 

MV 

tion  €tOc,t|r)  of  the  pressure  field  to  cKx.tlr)  is  the  negative  of  the 
pressure  gradient  in  the  fluid  element  at  time  r.  Except  in  degenerate 
cases , 

o(x,t|r)  f  -J7p(x,t|r)  (t  i  r),  (8.39) 

because  x  is  a  Lagrangian  coordinate  at  time  r,  not  a  laboratory  Cartesian 
coordinate.  Consequently,  the  Lagrangian  spacetime  structure  function 
^[p(x,t)  -  p(x' ,t|r)]2^  is  not  obtained  by  simply  extending  ( 8. 33)— ( 8 . 38) 
to  r  i  t. 


9.  COMPARISON  OF  INERTIAL-RANGE  PREDICTIONS  WITH  EXPERIMENT 

The  qualitative  inertial-range  predictions  of  the  abridged  LHDI 
equations  are  implied  by  Kolmogorov's  hypotheses  and  the  associated  dimen¬ 
sional  analysis.  Experimental  testing  of  the  results  therefore  involves, 
first  of  all,  a  test  of  the  Kolmogorov  theory.  The  Kolmogorov  scaling 
laws  have  received  support  from  measurements  of  energy  spectra  at  high 
Reynolds  numbers,  and  such  measurements  appear  to  offer  the  principal 
check  available  on  quantitative  predictions  of  inertial-  and  dissipation- 
range  properties. 

The  Kolmogorov  constant  is  customarily  estimated  from  experiment  by 
choosing,  visually,  best-fitting  straight  lines  to  log-log  plots  of  spectrum 
versus  wavenumber.  The  normalization  parameter  e  is  determined  by  drawing 
a  best-looking  curve  through  a  plot  of  the  dissipation  spectrum,  and  using 
isotropic  relations.  The  fitting  of  a  best  straight  line  requires  an  ex¬ 
plicit  or  implicit  estimate  of  the  maximum  wavenumber  at  which  dissipation 
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effects  on  spectrum  shape  can  be  neglected  in  comparison  with  experimental 
scatter.  This  estimate  plays  a  significant  role  because,  as  the  Figures  to 
follow  illustrate,  the  experiments  cover  very  limited  inertial  ranges. 
Consequently,  it  is  desirable  to  compare  the  data  directly  with  predictions 
of  the  entire . Kolmogorov  spectrum,  instead  of  merely  comparing  the  experi¬ 
mental  and  theoretical  estimates  of  the  Kolmogorov  constant. 

The  abridged  LHDI  prediction  of  the  Kolmogorov  spectrum  in  the  dissipa¬ 
tion  range  has  been  obtained  by  taking  U(k |-tQ )  oe  over  a  wide  k 

range,  and  then  integrating  (2.6)-(2.8)  forward  in  time,  with  nonzero  vis¬ 
cosity.  The  calculation  was  continued  until  the  normalized  dissipation 
2  1/4  5/4 

spectrum  (kAg)  E(k)/(e  v  )  reached  and  held  steady-state  values. 

Further  details  of  the  calculation  are  given  in  the  Appendix.  One-dimen¬ 
sional  spectra  were  then  computed  from  the  formulas 

•»  00 

^(k)  =  i  /  E(p)(l  -  k2/?2)?'^?,  *tot(k)  =  /  E(p)p"’1dp.  (9.1) 

k  k 

12 

Here  ♦^(k)  is  the  longitudinal  one-dimensional  spectrum  and  $  (k)  is 

the  sum  of  the  longitudinal  and  the  two  transverse  spectrum  functions.  The 
corresponding  asymptotic  inertial-range  spectra  are 

^(k)  *  (9/55)Ce2/3k"5/3f  *tot(k)  =  <3/5)Cc2/3k"5/3.  (9.2) 

The  computed  values  of  (kA  J^fjOOAe^v5^4)  from  the  inertial  and 
dissipation-range  calculations  are  plotted  together  in  Fig. 10.  This  kind 
of  plot  displays  differences  in  inertial-range  spectra  more  prominently  than 
the  usual  log-log  spectrum  plot.  Note  that  the  low-wavenumber  end  of  the 
decay-calculation  curve  rises  slightly  above  the  asymptotic  inertial-range 
line  C  =  1.77.  This  is  due  to  truncation  of  the  wavenumber  range  and  to 
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the  finite  time  of  evolution,  each  of  which  tends  to  depress  e  below  its 

asymptotic  value.  A  more  precise  computation  would  yield  a  curve  lying 

slightly  below  that  shown,  in  the  range  .002  <  (k/kg)  <  .02. 

The  experimental  points  plotted  in  Fig.  10  are  sea-water  data  obtained 

by  Grant,  Stewart,  and  Moilliet.  These  experiments  appear  to  be  the 

highest-Reynolds-number  measurements  of  absolute  spectrum  levels  which 

have  been  reported  (R^  <v.  3000).  The  points  represent  the  three  runs  on 

3  October  1959.  The  run  at  0950  hrs  was  cited  as  particularly  satisfactory 

with  regard  to  scatter  and  noise  levelj  it  is  plotted  in  different  ways  in 

o  1/4  5/4 

Figs.  8-10  of  Ref.  18.  The  dissipation  spectra  (k/kg)  *1(k)/{ e  v  ) 

for  the  three  runs  are  displayed  in  a  linear  plot  in  Fig,  11,  together 

with  the  theoretical  curve.  Both  data  and  theory  indicate  a  maximum  at 

about  (k A  )  *  0.09.  At  this  value  of  k/k  ,  the  curve  in  Fig.  10  is  down 
s  a 

from  the  asymptotic  inertial-range  line  by  about  30%. 

Figure  10  illustrates  the  difficulties  in  accurate  experimental 
determination  of  C.  Only  eight  points  lie  in  the  range  (k/kg)  <  .01, 
where  the  theoretical  curve  indicates  negligible  deviation  from  the  asymp¬ 
totic  inertial-range  level.  These  eight  points  scatter  from  20%  below  to 
25%  above  the  computed  inertial-range  line.  The  17  runs  reported  in  Ref. 

18  show  substantially  greater  scatter  at  (k/ks>  *  than  in  the  dissipa¬ 
tion  range  .01  <  (k/k  )  <  .5,  both  individually  and  collectively.  The 
fluctuations  at  low  k  seem  oscillatory  rather  than  random,  possibly  as  a 
result  of  coupling  between  the  turbulence  and  the  towed  body  on  which  the 
sensors  were  mounted.  Figure  12  compares  the  three  longest  runs  on  1 
February  1960,  and  on  2  February  1960,  with  the  theoretical  curve.  The 
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oscillations  are  particularly  prominent  in  the  run  at  1331  hrs  on  1  February. 

Grant,  Stewart,  and  Moilliet  report  a  mean  determination  C  =  1.44  from 

the  17  runs.  The  rms  fluctuation  of. the  17  values  of  C  about  the  mean  is 

about  0.19.  The  accuracy  of  absolute  level  determination  is  reported  as 

10%.  The  mean  value  of  C  rises  if  best  visual  fits  are  made  to  curves 

19 

like  that  in  Fig.  10  instead  of  simply  to  straight  lines.  This  is  be¬ 
cause  the  curve  approaches  the  asymptotic  inertial-range  line  gradually. 

3 

Heisenberg's  transfer  expression  gives  a  sharper  knee.  The  mean  value 
rises  further  if  the  longer  runs  are  given  greater  weight,  as  Figs.  10  and 
12  suggest . 

Probably  the  following  conclusions  are  justified:  (1)  The  experimental 

data  are  consistent  with  Kolmogorov's  theory.  (2)  The  data  suggest  that 

the  abridged  LHDI  prediction  for  C  is  better  than  an  order-of-magnitude 

approximation  to  the  correct  value.  (3)  The  abridged  LHDI  value  appears 

to  be  high  rather  than  low,  but  the  data  do  not  set  a  lower  bound  to  the 

magnitude  of  the  error  nor  determine  its  sign  unequivocally.  On  the 

theoretical  side,  the  approximations  made  in  deriving  the  LHDI  equations 

are  sufficiently  drastic  that  if  theory  and  experiment  agreed  to  better 

than  the  order  of  10%  it  would  be  an  accident.  An  interesting  question 

is  how  much  the  abridged  LHDI  value  of  C  differs  from  the  prediction  of 

the  unabridged  equations.  Numerical  integration  of  the  latter  in  the 

inertial  range  appears  to  be  feasible. 

Figure  13  shows  the  computed  values  of  (k/k  )5//3$  .  <k)/( 

8  TOT 

20 

plotted  together  with  data  obtained  by  M.  M.  Gibson  in  a  round  air  jet 
at  750.  The  data  were  normalized  by  values  of  e  determined  by  Gibson 
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from  the  values  only.  For  the  off-axis  run,  this  e  agrees,  to  within 

the  accuracy  of  determination,  with  the  e  determined  graphically  from  the 

4  .  values.  For  the  on-axis  run,  the  $  values  give  a  larger  e,  but 

tot  toi 

they  yield  so  irregular  a  dissipation  spectrum  that  accurate  determination 

is  impossible.  Therefore  the  on-axis  points  for  ♦  should  be  displaced 

19 

downward,  along  lines  of  slope  -4/3,  by  an  uncertain  amount. 

Gibson  estimates  C  *  1.57  for  the  on-axis  run  and  C  =  1.62  for  the 
off-axis  run,  on  the  basis  of  visual  straight-line  fits  to  the  <t>  ^  measure¬ 
ments  only.  As  in  the  sea-water  experiments,  these  values  rise  if  fits 
are  made  to  scaled  theoretical  curves  instead.  Both  runs  exhibit  aniso¬ 
tropy  in  the  inertial  range.  For  this  reason,  and  because  represents 
only  3/11  of  the  energy  in  an  isotropic  -5/3  range,  it  seems  more  justified 
to  compare  with  isotropic  theory,  rather  than  But  the  argument 
loses  force  if  the  transverse  spectrum  measurements  are  less  reliable 
than  the  longitudinal. 

It  should  be  stressed  in  conclusion  that  no  reported  experiments 
-5  /3 

confirm  the  k  asymptotic  spectrum  law  beyond  reasonable  doubt.  Actually, 

the  data  at  k/k  <  .01  on  Figs.  10  and  12  fit,  say,  k'7^4  better  than  k-5^3, 

3 

-5 /3 

but  the  scatter  is  so  great  that  no  inconsistency  with  k  can  be  inferred. 

Clearly  an  inertial  range  of  many  decades  extent  is  needed  for  a  definitive 

—5 /3 

answer.  Corrections  of  logarithmic  type  to  the  k  law  would  be  parti¬ 
cularly  difficult  to  resolve. 
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APPENDIX.  NUMERICAL  METHODS 

The  numerical  integrations  are  based  on  a  scheme  described  in  detail 

Q 

previously.  The  wavenumber  integrations  in  (2.6)-(2,8)  are  performed 
just  as  in  Ref.  9.  For  decay  studies,  (2.6)-(2.8)  are  integrated  forward 
in  time  from  tQ  by  a  modification  of  the  time-integration  scheme  of  Ref.  9: 
The  first  two  terms  on  the  right-hand  sides  of  (2.7)  and  (2.8)  are  treated 
in  parallel  with  the  vk  terms  on  the  left-hand  sides.  Thus, for  (2.7),  a 
function  N(k;t|r)  is  defined  by 

t 

N(k;t|r)U(k;t|r)  =  //ACCkp<iU(k;t|r)-DkpqU(p;t|r)]/  U(q;t|s)ds,  (Al) 

2 

and  the  contribution  of  the  two  terms  is  included  by  replacing  vk  with 
vk2  +  N(k;t|r)  in  the  scheme  of  Ref.  9.  N(k;t|r)  is  evaluated  at  the 
mid-point  of  each  time  step,  like  the  other  quantities  in  the  scheme. 

The  analogous  procedure  is  used  for  (2.8).  This  modification  gives  improved 
stability  and  accuracy.  The  C  and  D  terms  on  the  right-hand  side  of  (Al)  are 
treated  in  strict  parallel  with  each  other  to  avoid  errors  from  the  effect 
of  time  and  wavenumber  discretizing  upon  the  cancellations  of  the  small-q 
contributions . 

For  the  steady-state  solution  of  (5.4)  and  (5.5),  this  method  is 

augmented  by  an  iteration  procedure.  The  time  limit  “  in  (5.5)  is  truncated 

_o 

to  a  large  finite  value  MAt  (M  integer).  Initial  guesses  G  (t)  and  R  (t) 
are  taken  (identical  Gaussian  functions).  The  values  (f(0)  =  1,  R(0)  =  1 
are  held  fixed  and  improved  functions  obtained  by  a  modification  of  the 


scheme 
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G^t  )  =  -  {Gr“1(T  )  +  g’V  .)  +  i  AtClTC t  )  +  H(t  )]},  !T(0)  =  0, 
n  2  n  n-i  2  n  n-x 

R^t  )  =  4  {Rr’1(T  )  +  .)  +  i  At[S(i  )  +  S(t  .)]},  2‘(0)=0.  (A2) 

n  i  n  n-i  z  n  n-x 


Here  G* (Tn>*  Fr(Tn)  rth  approximation  at  the  nth  time  step*  and 

H(t),  S"(t)  are  the  right-hand  sides  of  (S.4),  (5.5).  For  each  r,  (A2)  is 

used  for  each  n  in  succession,  starting  with  n  =  1  [t^  =  nAt  (n  =  1,2, 

time—  _  __ 

...M)].  At  each  step,  S  and  H  are  evaluated  by  the  trapezoidal  rule  of  Ref. 


9,  using  always  the  best  available  values  of  G’  and  Rs  As  soon  as  G  (in) 

and  R1* ( xn)  are  computed  for  any  r  and  n,  the  old  values  Gr  1(tn)  and  R*1  t^) 

_  2/3 

are  discarded  completely.  Functions  like  R(w  s)  are  evaluated  from  the 
R(tn)  by  two-point  interpolation.  The  modification  of  (A2),  in  accord 
with  (Al)  seq.,  is  to  replace  ^Cs’(rn)  +  Tn_1> 3  as  follows: 

[l  -  exp{-  •|[N(xn)  +  N(in_1)]At} 


7  At[S(x  )  +  S(T  )] 
i  n  n-x 


N(t  )  +  N(t  .) 
n  n-l 


X  [S(t  )  +  S(t  .)],  (A3) 

n  n-x 

where  -IT(t)R(t)  is  defined  as  the  first  two  terms  on  the  right-hand  side  of 
(5.5).  The  (T  equation  is  treated  in  parallel  fashion. 

The  iteration  scheme  was  found  to  converge  reliably  but  slowly;  e.g., 
Ir1^)  -  Rr_1(Tn)|  decreased  typically  by  20%  per  itjaration.  Convergence 
is  markedly  speeded  by  repeated  cycles  of  iteration-followed-by-extrapolation. 
The  extrapolation  assumes  geometrically  decreasing  errors  (the  errors  do 
decrease  nearly  geometricallly). 


-64- 


Equation  (5,11)  is  integrated  by  the  same  wavenumber  discretizing 
used  in  (5,4)  and  (5,5).  The  contributions  from  wavenumber  triads  with 
different  a  values  are  accumulated  separately  so  as  to  yield  W(a)  without 
the  need  for  further  calculation. 

The  numerical  values  reported  in  the  text  come  from  a  calculation  in 
which  the  v  and  w  integrations  in  (5.4),  (5.5)  extend  over  the  range 
2*145/12  ^  2^3  x  io”1*)  to  2^®^  ('v  96.),  which  is  divided  into  112  one- 
sixth-octave  steps  according  to  the  method  of  Ref.  9.  The  contributions  of 
higher  w  and  v  values  are  approximated  by  correction  terms  based  on  (4.13); 
these  corrections  are  quite  small.  The  s  integrations  are  truncated  at 
s  -  20,  and  R(t),  G(t)  are  computed  over  (0  <  t  £  20)  at  steps  of  0.1. 
Iteration  was  continued  until  the  change  in  C  per  iteration  was  <0.1%,  with 
extrapolation  indicating  residual  error  <0.3%.  Changes  of  R(s)  and  G(s)  per 
iteration  were  less  than  both  .001  and  1%,  with  extrapolation  indicating 
residual  errors  less  than  .005  and  5%.  Errors  from  finite  step-sizes  and 
truncation  were  estimated  by  varying  the  relevant  parameters  and  extrapolating. 
The  estimated  total  error  in  C  from  all  causes  is  about  1%.  Errors  in  some 
of  the  other  dimensionless  constants  are  larger,  but  are  not  expected  to 
exceed  3%.  The  theory  probably  does  not  justify  more  accurate  computation. 

Six  to  eight  decimal  digits  were  carried  in  the  computations,  and  all 
final  numbers  reported  in  the  text  are  rounded  to  two  or  three  digits.  To 
six  digits,  the  values  found  for  (6.1)  and  (6,2)  are  C  =  1.76612, 

IQ  s  1.04404,  IR  =  1.87438,  1^  =  .760973,  =  -.188401.  The  values  com¬ 

puted  for  the’  left  and  right  sides  of  (8.26)  are  8.73360  and  8.96029,  which 
agree  to  within  2.6%, 
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For  rough  calculations,  the  computed  functions  G(t)  and  R(x)  are 
adequately  approximated  by 

G(t)  =  4  exp(2T/IG)/[l  +  exp(2t/IG)]2,  R(t)  =  [l+(j  (A4) 

The  decay  calculation  used  in  Sec.  9  was  performed  by  taking  v  =  .008, 
E(k,tQ)  =  6.29  k“5^3,  and  a  k  range  2-1976  .11)  to  23776  72.)  dis¬ 

cretized  into  28  one-third  octave  steps.  Equations  (2.6)-(2.8)  were  in¬ 
tegrated,  with  the  techniques  described  above,  for  20  time  steps  At  =  .02, 
by  which  time  changes  in  the  normalized  energy  and  dissipation  spectra 
were  imperceptible  on  automatically  generated  plots.  Final  values  kg  =  58.5 
and  =  440  were  obtained,  with  respective  changes  of  0.1%  and  0.2%  per 
time  step.  The  estimated  overall  computational  errors  are  not  significant 
for  the  comparison  with  experiment  in  Figs.  10-13,  provided  k/kg  <0.6. 

The  difference  between  this  Rx  =  440  spectrum  and  an  =  ■  spectrum 

normalized  in  the  same  way  probably  is  negligible  also.  It  was  to  minimize 

-5 /3 

the  difference  that  the  initial  spectrum  was  taken. <*  k  down  to  the  low¬ 
est  wavenumber  retained. 

All  the  calculations  were  performed  on  the  IBM  7094  computer  at  the 
Goddard  Institute  for  Space  Studies.  The  plots  were  generated  on  a  General 
Dynamics  S-C  4020  Recorder. 
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Peterborough,  New  Hampshire 

The  minus  signs  should  be  changed  to  plus  signs  in  the  third  line  of 

(8.10),  the  first  line  of  (8.11),  and  the  second,  fifth,  and  seventh  lines 

of  (8.12).  In  the  sixth  line  of  (8.11),  p  should  be  c  .  In  the  fourth 

c  •  c 

line  of  (10.9),  the  plus  sign  should  be  a  minus  sign. 


ERRATA:  PRELIMINARY  CALCULATION  OF  THE  KOLMOGOROV  TURBULENCE  SPECTRUM 
[Phys.  Fluids  8,  995  (1965)] 

Robert  H.  Kraichnan 
Peterborough,  New  Hampshire 

In  the  third  line  of  (3),  the  plus  sign  should  be  a  minus  sign.  This 
error  was  carried  through  the  numerical  computations.  The  corrected  numer¬ 
ical  coefficients  in. (10)  are  1.60,  0.50,  0.99,  in  order  of  appearance.  The 
corrected  coefficients  in  (10*)  are  1,77,  0.53,  0.96.  Corrections  in  the 
Figures  are  small.  For  corrected  Figures,  see  R.  H.  Kraichnan,  Phys. 

Fluids  9,  xxx  (1966)  [Research  Report  8], 
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FIGURE  CAPTIONS 

Fig.  1.  Integration  paths  for  following  the  evolution  of  the  Lagrangian 
velocity. 

Fig.  2.  The  inertial-range  Lagrangian  time-correlation  and  response 
functions  R(t)  and  G(t),  versus  nondimensional  time  difference  t. 

Fig.  3.  Dimensionless  measures  of  localness  of  inertial-range  energy 
transfer, versus  maximum  to  minimum  wavenumber  ratio  a.  Here  W(a)  =  W(a)/e 
and  P(a)  =  P(a)/e. 

Fig.  4.  Dimensionless  Lagrangian  spacetime,  or  two-particle,  structure 
functions  F(s)  and  F(|  (s)  in  the  inertial  range,  versus  dimensionless  time 
difference  s.  The  straight  line  is  the  function  ics  of  (7.12a). 

Fig.  5.  Dimensionless  Eulerian  spacetime  structure  functions  FE(s)  and 
Fj^(s)  in  the  inertial  range,  versus  dimensionless  time  difference  s. 

Fig.  6.  Nondimensional  measures  F(s)  and  F()  (s)  of  the  inertial-range 
covariance  of  one  particle's  velocity  with  a  second  particle's  prior 
acceleration,  versus  dimensionless  time  difference  s. 

Fig.  7.  Nondimensional  time-displaced,  inertial-range  cospectrum  Q(s)  of 
pressure  gradient  and  prior  Lagrangian  velocity,  versus  nondimensional 
time  difference  s. 

Fig.  8.  Nondimensional  inertial-range  covariances  Z(s)  and  Zx(s)  of  one 
particle's  prior  velocity  with  a  second  particle's  acceleration,  versus 
dimensionless  time  difference  s. 

Fig.  9.  Nondimensional  inertial-range, two-particle, Lagrangian  acceleration 
covariances  Z(s)  and  Zx(s),  versus  dimensionless  time  difference  s. 
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Fig.  10.  Computed  nondimensionalized  spectrum  function  k  <(;  ( k ) 

in  the  inertial  and  dissipation  ranges  compared  with  the  October  1959  data 

of  Grant,  Stewart,  and  Moilliet.  Circles,  run  0905/3/10/59;  dots, 

0907/3/10/59;  plus  signs,  0915/3/10/59.  (The  low-wavenumber  terminus  of 

the  decay  curve  lies  within  the  fourth  circle  from  the  left.) 

Fig.  11.  Computed  nondimensionalized.  dissipation  spectrum  k  4>^< k >  compared 

with  the  October  data  of  Grant,  Stewart,  and  Moilliet.  Circles,  run 

0905/3/10/59;  dots,  0907/3/10/59;  plus  signs,  0915/3/10/59, 

5 /3 

Fig.  12(a).  Computed  function  k  ^(k)  compared  with  the  three  longest 

runs  obtained  by  Grant,  Stewart,  and  Moilliet  on  1  February  1960.  Circles, 

run  1126/1/2/60;  dots,  1203/1/2/60;  plus  signs,  1331/1/2/60. 

5 /3 

Fig.  12(b),  Computed  function  k  <f^(k)  compared  with  the  three  longest 

runs  obtained  by  Grant,  Stewart,  and  Moilliet  on  2  February  1960,  Circles, 

run  1301/2/2/60;  dots,  1316/2/2/60;  plus  signs,  1445/2/2/60. 

Fig.  13.  Comparison  of  computed  spectra  with  data  of  M.  M.  Gibson.  Upper 

5/3 

curve  and  data:  computed  function  k  ^tot^’  <*ots»  on-axis  run;  plus 

5 /3 

signs,  off-axis  run.  Lower  curve  and  data:  computed  function  k  4>1( k ) ; 

dots,  on-axis  run;  plus  signs,  off-axis  run.  The  data  for  ♦  appear  only 
at  those  wavenumbers  where  Gibson  obtained  the  necessary  measurements  of 
both  longitudinal  and  transverse  spectra. 
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